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PREFACE 

The  work  described  in  this  report  was  authorized  by  Headquarters,  US  Army 
Corps  of  Engineers  (HQUSACE) ,  as  part  of  the  Geotechnical-Rock  Problem  Area 
of  the  Repair,  Evaluation,  Maintenance,  and  Rehabilitation  (REMR)  Research 
Program.  The  work  was  performed  under  Work  Unit  32306,  "Stability  of  Large 
Concrete  Structures  on  Rock,"  for  which  Mr.  James  B.  Warriner  was  Principal 
Investigator.  Mr.  Lewis  Gustafson  (CECW-EG)  was  the  REMR  Technical  Monitor 
for  this  work. 

Mr.  Jesse  A.  Pfeiffer,  Jr.  (CERD-C)  was  the  REMR  Coordinator  at  the 
Directorate  of  Research  and  Development,  HQUSACE;  Mr.  James  E.  Crews  (CECW-OM) 
and  Dr.  Tony  C.  Liu  (CECW-ED)  served  as  the  REMR  Overview  Committee;  Mr. 
William  F.  McCleese  (CEWES-SC-A) ,  US  Army  Engineer  Waterways  Experiment  Station 
(WES)  was  the  REMR  Program  Manager.  Mr.  Jerry  S.  Huie  was  the  Problem  Area 
Leader. 

The  work  was  performed  at  the  University  of  Missouri-Rolla  and  this  report 
was  prepared  by  Dr.  James  R.  Carr,  Department  of  Geological  Sciences  under 
an  Intergovernmental  Personnel  Act  Agreement  and  under  the  general  supervision 
of  Dr.  William  F.  Marcuson  III,  Chief,  Geotechnical  Laboratory  and  the  direct 
supervision  of  Dr.  Don  C.  Banks,  Chief,  Soil  and  Rock  Mechanics  Division. 

Commander  and  Director  of  WES  was  COL  Larry  B.  Fulton,  EN.  Technical 
Director  was  Dr.  Robert  W.  Whalin. 
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CONVERSION  FACTORS,  NON-SI  TO  SI  (METRIC) 
UNITS  OF  MEASUREMENT 


Non-SI  units  of  measurement  used  in  this  report  can  be  converted  to  SI 
(metric)  unit  as  follows: 

Multiply  By  To  Obtain 

feet  0.3048  metres 
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SURFACE  ROUGHNESS  CHARACTERIZATION  OF  ROCK  MASSES 
USING  THE  FRACTAL  DIMENSION  AND  THE  VARIOGRAM 


PART  I:  INTRODUCTION 

1,  This  report  summarizes  the  first  known  attempt  to  describe  the 
roughness  of  rock  surfaces  using  the  concepts  of  the  fractal  dimension  and  the 
variogram.  For  this  investigation,  joint  surfaces  adjacent  to  Libby  Dam,  MT, 
were  used.  One  joint  surface  was  part  of  a  rock  slide  designated  as  the 
DS  +  122  slide.  The  objective  of  this  study  was  to  examine  the  utility  of  the 
fractal  dimension  and  the  variogram  for  the  classification  of  surface  rough¬ 
ness  of  rock  masses. 


General  Description  of  the  DS  +  122  Rock  Slide 


2.  The  DS  +  122  rock  slide  occurred  at  Libby  Dam,  MT,  on  January  31, 

1971  (Hamel  1974)  .  This  slide  consisted  of  a  rock  wedge  which  slid  toward  the 

3 

dam  abutment.  The  volume  of  this  rock  wedge  was  approximately  46,000  m  .  The 
DS  +  122  wedge  was  formed  and  bounded  by  the  intersection  of  a  bedding  plane 
and  a  fracture  (joint).  From  a  qualitative  viewpoint,  the  surface  of  the 
bedding  plane  was  quite  smooth  whereas  the  joint  surface  had  a  rougher  texture 
(Banks  and  Strohm  1974).  This  report  presents  a  quantitative  description  of 
the  roughness  of  the  joint  surface. 

3.  This  description  involved  the  characterization  of  the  spatial  varia¬ 
tion  in  the  roughness  across  the  planar  surface  of  this  joint.  In  analyses  of 
the  sliding  stability  of  rock  masses,  asperities  along  a  potential  slip  sur¬ 
face  are  a  key  consideration  (Patton  1966;  Patton  and  Deere  1970;  Hoek  and 
Bray  1974;  and  Goodman  1976).  This  report  does  not  describe  the  characteris¬ 
tics  of  the  asperities  along  the  DS  +  122  joint  surface.  A  quantitative 
description  of  the  roughness  of  this  slide  surface  derived  from  a  field  survey 
is  presented.  From  this  description,  it  is  hoped  that  a  method  for  deriving 
properties  of  slip  surface  asperities  can  be  developed. 
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Fractal  Dimensions  and  Variograms 


4.  A  quantitative  assessment  of  the  roughness  of  the  DS  +  122  slide 
surface  was  completed  using  the  concept  of  the  fractal  dimension.  This  con¬ 
cept  describes  the  degree  of  variation  a  natural  surface  has  from  its  topolog¬ 
ical  facsimile.  Topology  is  a  branch  of  mathematics  devoted  to  the  study  of 
geometrical  configurations.  For  example,  the  topological  model  of  a  joint  in 
a  rock  is  a  plane.  In  reality,  a  joint  is  not  a  perfect  plane,  but  has  sur¬ 
face  irregularities  which  impart  a  greater  surface  area  to  the  joint  when 
compared  to  its  topological  model.  The  fractal  dimension  describes  the  degree 
of  variation  a  curve,  a  surface,  or  a  volume  has  from  its  topological  ideal. 

5.  For  example,  a  coastline  is  modeled  topologically  as  a  line  or  a 
semi-circle.  The  topological  dimension  is  therefore  one  (a  one  dimensional 
feature).  A  coastline  is  rarely  linear  or  circular,  however,  and  is  often 
highly  irregular.  This  irregularity  gives  the  coastline  a  greater  length 
compared  with  its  topological  model.  Hence,  the  fractal  dimension  of  such  an 
irregular  coastline  is  greater  than  one.  A  fractal  is  a  curve,  a  surface,  or 
a  volume  whose  dimension  exceeds  its  topological  dimension.  The  more  irregu¬ 
lar  or  rough  a  surface  is,  the  greater  will  be  its  fractal  dimension. 

6.  Once  the  irregularity  of  a  surface  has  been  described  by  the  fractal 
dimension,  the  spatial  variation  of  the  surface  irregularities,  or  asperities, 
can  be  described  using  variograms.  The  variogram  describes  the  spatial  con¬ 
tinuity  and  relationship  inherent  to  a  type  of  regionally  distributed  natural 
data.  The  variogram  can  be  used  to  show  the  spatial  behavior  of  the  asperi¬ 
ties  whether  these  irregularities  are  of  a  local  extent  with  little  spatial 
continuity  or  are  continuous  with  a  more  regular  spatial  variation.  Further¬ 
more,  the  variogram  can  be  used  to  evaluate  the  different  magnitudes  of  asper¬ 
ity  sizes  which  are  present  on  a  rock  surface.  The  resolution  of  these 
asperities  is  important.  Goodman  (1976)  has  shown,  for  example,  that  asperi¬ 
ties  increase  the  friction  angle  of  a  joint  by  an  amount  equal  to  the  mean  or 
predominate  asperity  angle.  Hence,  any  investigation  of  the  sliding  stability 
of  a  joint  (or,  in  general,  any  rock  surface)  must  be  concerned  with  the 
characteristics  of  the  asperities  along  the  surface. 


9 


PART  II:  FRACTAL  DIMENSION  OF  JOINT  SURFACES 


Fractal  Dimension 

7.  A  detailed  description  of  fractals  and  the  fractal  dimension  is  pro¬ 
vided  in  Appendix  A.  In  summary,  a  fractal  is  a  line,  a  surface,  or  a  volume 
whose  Hausdorf f-Besicovitch  dimension  (i.e.,  fractal  dimension,  see  Appen¬ 
dix  A)  exceeds  its  topological  dimension.  In  Appendix  A,  the  coastline  of 
Great  Britain  is  used  as  an  example.  The  topological  dimension  of  this  coast 
is  1  (a  line  or  a  semi-circle) .  The  length  of  this  coast  can  be  measured  by 
selecting  a  particular  segment  length,  y  ,  then  counting  the  number,  N  ,  of 
these  y-length  segments  which  are  needed  to  approximate  the  total  length,  L  , 
of  the  coastline.  This  length  is  computed  as 

L  =  Ny  (1) 

But,  as  we  are  realizing  when  we  measure  with  increasing  resolution,  we 
encounter  an  increasing  number  of  microvariabilities,  thus  equation  (1)  should 
be  expressed  as 

L  =  Ny  (la) 

This  is  a  phenomenon  that  has  no  real  mathematical  proof  and  is  conceptually 
shown  in  Appendix  A. 

8.  It  is  an  interesting  phenomenon  that  as  y  decreases,  L  increases 
without  bound  (for  the  actual  coast,  not  for  its  mapped  representation;  a  map 
has  a  certain  resolution  limit  and  increments  of  y  below  that  limit  yield 
the  same  measure  of  length,  L  ,  hence  with  a  map,  L  does  have  a  limit). 
Despite  this  fact,  if  we  plot  log^(N)  versus  log^(y),  a  straight  line  is 
obtained.  The  negative  of  the  ~lope  of  this  line  is  the  Hausdorf f-Besicovitch 
dimension,  D  .  This  dimension  is  the  fractal  dimension  (Mandelbrot  1982). 

9.  The  fractal  dimension  can  be  incorporated  into  equation  (1)  to 

become 

L  =  Ny°  (2) 


10 


If,  throughout  the  plot  of  logj^(N)  versus  log^Cy),  D  is  found  to  be  greater 
than  the  topological  dimension,  a  fractal  is  defined.  In  Appendix  A,  for 
instance,  the  fractal  dimension  of  the  west  coast  of  Great  Britain  is  shown  to 

be  1.25.  This  is  greater  than  the  topological  dimension,  1  ;  hence,  the  west 

coast  of  Great  Britain  is  a  fractal. 

10.  For  the  analysis  of  the  surface  roughness  of  the  joint  bounding  the 

DS  +  122  slide,  the  concept  of  the  fractal  dimension  is  intuitively  useful. 

The  more  irregular  a  surface  is,  the  larger  will  be  its  fractal  dimension. 

The  DS  +  122  joint  surface,  however,  is  planar  with  a  topological  dimension  of 
two  (D^  =2).  To  measure  the  fractal  dimension  of  this  surface,  equation  (2) 
is  rewritten  to  become 


A  =  Na°  (3) 

where  A  is  the  total  area  of  the  surface  and  a  is  some  small,  regular  area 
analogous  to  y  in  equation  (2).  From  equation  (3),  D  is  found  by  plotting 
log^(N)  versus  log^(a).  For  two  dimensional  surfaces,  a  fractal  surface  is 
defined  if  D  exceeds  2. 

Cross-Section  Segmentation  of  the  DS  +  122  Joint  Surface 

11.  For  the  elevation  of  the  DS  +  122  slide  surface,  the  solution  for 
D  using  equation  (3)  is  difficult  to  obtain.  As  an  alternative,  cross  sec¬ 
tions  were  plotted  along  discrete  sections  of  the  joint  surface.  These  cross 
sections  were  plotted  using  a  contour  plot  derived  from  ground-based  photo- 
grammetry.  The  map  scale  was  1:40  and  the  contour  interval  was  2  cm.  In 
total,  thirty  cross  sections  were  analyzed.  Of  these,  twenty  were  orthogonal 
to  the  direction  of  shear  failure,  and  ten  were  parallel  to  the  direction  of 
shear  failure.  The  locations  of  these  thirty  cross  sections  are  shown  in 
Figure  1. 

12.  Selected  cross  sections  are  plotted  in  Figures  2  and  3.  By  analyz¬ 
ing  cross  sections,  the  two-dimensional  analysis  offered  by  equation  (3)  was 
transformed  to  thirty  one-dimensional  analyses  using  equation  (20).  Topologi¬ 
cally,  the  trace  of  each  cross  section  is  a  line.  A  particular  cross  section 
was  determined  to  be  a  fractal  if  its  fractal  dimension  was  greater  than  1. 
Once  calculated  for  each  cross  section,  the  fractal  dimension  served  as  an 


11 


DS  4-  122  joint  surface  showing  cross-section  locations 


indicator  of  surface  roughness.  Because  each  cross  section  represented  a  uni¬ 
que  spatial  location  along  the  DS  +  122  slide,  the  analysis  of  thirty  cross 
sections  described  the  spatial  variation  in  the  surface  roughness  along  the 
slide. 

Fractal  Dimension  of  Selected  Cross  Sections 

13.  A  fractal  dimension  was  computed  for  each  of  the  thirty  cross  sec¬ 
tions  plotted  in  Figure  1.  To  minimize  the  length  of  this  report,  Table  1 
presents  calculations  only  for  those  cross  sections  shown  in  Figures  2  and  3. 
Calculations  of  fractal  dimensions  for  other  cross  sections  and  presented  in 
Appendix  C. 

Table  1 

Fractal  Dimension  Calculations  (Method  2) 


Cross 

Number  of 

Segment  Increments 

,  N 

Section 

y  =  305  cm 

y  =  100  cm  . 

y  =  20  cm 

-slope  =  D* 

AD 

4.056 

12.400 

62.500 

1.0038 

AE 

3.810 

11.680 

58.900 

1.0050 

AL 

3.900 

9.240 

60.968 

1.0091 

AS' 

4.118 

12.610 

63.700 

1.0052 

BA 

7.793 

23.770 

118.843 

1.0000 

BF 

6.832 

20.940 

105.330 

1.0040 

*  Slope  =  (log10(N20)  -  log10<N305))/(log10(20)  -  log10(305)). 


14.  For  each  cross  section,  several  consistently  smaller  segment 
lengths,  y  ,  were  selected  for  the  calculation  of  the  fractal  dimension. 
Usually,  three  different  lengths  were  chosen:  305  cm,  100  cm,  and  20  cm. 

This  was  an  arbitrary  choice  and  seemed  to  work  well  for  this  application.  In 
no  way  do  these  lengths  conform  to  a  special  property  of  these  cross  sections. 
For  the  calculation  of  the  fractal  dimension,  measurable  lengths,  y  ,  must  be 
selected  by  arbitrary  decision.  For  this  study,  the  lengths  200  cm,  50  cm, 
and  5  cm  could  have  been  chosen  and  would  have  yielded  similar  results. 

15.  For  each  cross  section,  the  fractal  dimension,  D  ,  was  evaluated 
quantitatively  using  a  pair  of  segment  lengths,  y  .  In  addition,  log^CN) 
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versus  log^Cy)  plots  were  developed  for  each  cross  section  for  graphical 
analysis.  These  plots  are  shown  in  Figures  4  through  9  with  visually  best  fit 
lines.  The  negative  of  the  slope  of  each  plot  is  the  fractal  dimension.  The 
fractal  dimension  for  each  cross  section  was  determined  by  comparing  the 
logjg(N)  versus  log^(y)  graphical  results  to  the  quantitative  calculation 
results . 


Results  of  the  Fractal  Dimension  Calculations 


16.  The  results  of  the  fractal  dimension  calculations  for  all  thirty 
cross  sections  are  presented  in  Table  2.  To  make  these  results  more  meaning¬ 
ful,  Figure  10  shows  the  fractal  dimension  for  each  cross  section  plotted  on 
the  map  of  the  DS  +  122  joint  surface.  This  figure  demonstrates  the  spatial 
variation  in  the  fractal  dimension  across  the  DS  +  122  joint  surface. 

17.  Figure  10  shows  the  DS  +  122  joint  surface  between  cross  sec¬ 
tions  AB  through  AF  to  be  rather  smooth.  The  joint  surface  becomes  rougher 
toward  the  AZ  cross  section.  The  direction  parallel  to  sliding  has  a  rough¬ 
ness  which  increases  from  cross  sections  BA  and  B^  to  cross  sections  BI 

and  BJ.  The  direction  of  shear  failure  along  this  surface  was  from  cross 
sections  AZ,  BI ,  BJ  toward  cross  sections  AB,  BA,  BB. 

18.  It  is  Important  to  evaluate  the  method  used  to  calculate  the 
fractal  dimension  for  each  cross  section.  As  we  reviewed  previously,  the 
calculation  of  the  fractal  dimension  begins  by  selecting  a  segment  length, 

y  .  Then,  beginning  from  one  end  of  an  irregular  line,  the  number  of  steps, 

N  ,  is  counted  to  reach  the  opposite  end.  Here  a  dilemma  is  reached. 

In  counting  the  number  of  steps,  N  ,  what  happens  if  a  small  section  of  the 
irregular  line  remains  at  the  end  and  this  remainder  is  smaller  than  the 
length,  y  ?  It  is  necessary  to  add  the  remainder,  normalized  by  y  ,  to  the 
number  of  steps,  N  ;  i.e.,  N  =  N  +  remainder/y  (Mandelbrot  1985).  Table  2 
shows  that  this  remainder  is  considerably  important  for  calculations.  The 
method  1  calculation  does  not  consider  the  remainder  whereas  the  method  2 
calculation  does  include  the  remainder.  Hereafter,  in  this  study,  only  the 
method  2  calculation  is  used  to  determine  the  fractal  dimension  of  rock 
surfaces . 
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Figure 

4.  Log^CN)  versus 
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for  cross  section  AD 
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Figure  5.  Log^(N)  versus  log^Cy)  plot  for  cross  section  AE 
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Figure  6.  Log^(N)  versus  log^(y)  plot  for  cross  section  AL 
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Figure  7.  Log^(N)  versus  log^Q<y)  plot  for  cross  section  AS' 
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Figure  8.  Log^(N)  versus  log^(y)  plot  for  cross  section  BA 
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Figure  9.  Log^(N)  versus  log^Cy)  plot  for  cross  section  BF 
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Table  2 

Fractal  Dimension  Summary 


Cross  Fractal  Dimension 


Section 

Method  1* 

Method  2** 

AB 

1.020 

1.0027 

AD 

1.007 

1.0038 

AE 

1.060 

1.0050 

AF 

1.007 

1.0000 

AG 

1.020 

1.0072 

AH 

1.030 

1.0119 

AJ 

1.020 

1.0060 

AL 

1.070 

1.0100 

AN 

1.010 

1.0012 

AP 

1.010 

1.0009 

AQ 

1.030 

1.0061 

AR 

1.020 

1.0107 

AS 

1.030 

1.0100 

AS' 

1.030 

1.0052 

AT 

1.050 

1.0200 

AU 

1.010 

1.0040 

AW 

1.010 

1.0268 

AX 

1.050 

1.0189 

AY 

1.070 

1.0100 

AZ 

1.110 

1.0030 

BA 

1.040 

1.0000 

BB 

1 . 0u3 

1.0000 

BC 

1.007 

1.0000 

BD 

1.003 

1.0000 

BE 

1.040 

1.0000 

BF 

1.040 

1.0040 

BG 

:  030 

1.0000 

BH 

1.020 

1.0020 

BI 

1.010 

1.0030 

BJ 

1.040 

1.0050 

*  For  method  1,  only  whole  numbers  are  used  for  N 

**  For  method  2,  the  fractional  remainder  is  added  to  N 


Control  Calculation:  Fractal  Dimension  of  the 
Australian  Coast 


19.  A  control  calculation  was  sought  to  verify  the  fractal  dimension 
calculation  procedure.  In  a  previous  study,  Mandelbrot  (1967)  reviewed  the 
fractal  dimension  of  coastlines.  He  found  the  fractal  dimension  of  coastlines 
varied  from  a  maximum  of  1.25,  for  the  west  coast  of  Great  Britain,  to  a 
minimum  of  1.02,  for  the  South  African  coast.  For  the  Australian  coast, 
Mandelbrot  found  the  fractal  dimension  to  be  1.13.  This  work  is  a  standard 
reference  because  Mandelbrot  was  the  first  to  demonstrate  the  fractal  dimen¬ 
sion  of  coastlines;  moreover,  he  developed  the  concept  of  the  fractal 
dimension. 

20.  Mandelbrot's  results  give  interesting  examples  which  can  be  used 
for  control  calculations.  For  this  purpose,  the  coastline  of  Australia  was 
chosen.  The  outline  of  this  continent  is  shown  in  Figure  11.  Using  the  seg¬ 
ment  lengths,  y  ,  shown  in  Table  3,  the  fractal  dimension  was  found  to  be 
1.11  using  the  method  2  technique;  this  fractal  dimension  was  determined  from 
the  plot  shown  in  Figure  12.  This  result  is  close  to  the  value  reported  by 
Mandelbrot.  Mandelbrot  used  the  method  2  calculation  for  his  research 
(Mandelbrot  1985).  There  is  no  proof  that  has  been  presented  which  supports 
the  use  of  the  method  2  technique.  Its  use  is  justified  by  the  realization 
that  the  exclusion  of  the  remainder  can  introduce  considerable  error  in  the 
calculation  especially  if  the  segment  length,  y  ,  is  large  with  respect  to 
the  line  being  measured. 
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Figure  11.  Outline  of  the  coast  of  Australia 


Table  3 

Fractal  Dimension*  Calculation,  Australian  Coast 


Ruler 

length 

y  (km) 

lo810(y) 

Segment 

Counts,  N 

log10(N) 

500 

2.700 

24.442 

1.3881 

300 

2.477 

43.868 

1.6421 

100 

2.000 

154.471 

2.1889 

*  Fractal  dimension  is  1.11. 
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Figure  12.  Log^(N)  versus  log^Cy)  plot  for  the  Australian  coast 
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PART  III:  VARIOGRAM  ANALYSIS  OF  THE  DS  +  122  JOINT  SURFACE 


Variogram 


21.  In  brief  overview,  the  variogram  is  a  mathematical  function  which 
describes  the  spatial  similarity  between  the  value  of  a  geological  phenomenon 
at  one  location  and  the  value  at  a  different  location.  The  variogram  is 
described  in  detail  in  Appendix  B.  This  function  is  an  important  concept  in 
modern  geostatistics  and  has  the  form: 


2  Y(h)  =  1/N 


Z(x±)  - 


(A) 


where  Z  is  a  vector  containing  the  data  which  represents  the  spatial  pheno¬ 
menon  under  study,  N  is  the  number  of  pairs  of  these  data  separated  by  a 
distance,  h  ,  and  y  is  the  variogram  function.  By  describing  the  spatial 
similarity,  the  variogram  also  captures  spatial  continuity;  that  is,  the 
variogram  shows  how  rapidly  the  similarity  between  two  data  locations 
decreases  with  increasing  distance  separating  them. 

22.  For  rock  mechanics  applications,  the  variogram  function,  along  with 
geostatistics,  has  received  some  use.  These  applications  have  been  for  the 
analysis  of  the  spatial  variation  in  rock  mass  features,  particularly  fracture 
density  and  orientation  (LaPointe  1980).  Early  work  with  this  type  of  analy¬ 
sis  has  been  extended  to  the  simulation  of  rock  mass  properties  (Miller  and 
Borgman  1985).  The  variogram  has  also  been  used  to  describe  the  spatial  vari¬ 
ation  in  rock  mass  features  adjacent  to  and  involved  in  ground  subsidence  over 
abandoned  mines  (Van  Besien  1985).  Hence,  the  variogram  is  known  technique  in 
the  field  of  rock  mechanics. 


Variograms  of  Selected  Cross  Sections 


23.  Variogram  analysis  was  used  for  this  study  to  describe  the  spatial 
correlation  and  continuity  of  joint  surface  asperities.  The  objectives  of 
this  application  were  to:  (a)  simply  evaluate  the  utility  of  the  variogram 
function  for  this  type  of  analysis;  (b)  investigate  the  spatial  correlation 
and  continuity  of  rock  mass  asperities;  and  (c)  determine  the  correlation,  if 
any,  between  the  variogram  and  the  fractal  dimension. 
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24.  Selected  cross  sections  are  chosen  and  presented  for  t'  *s  applica¬ 
tion.  These  are  shown  in  Figures  13  and  14.  The  variograms  for  Lhese  cross 
sections  are  presented  in  Figures  15  through  20.  The  actual  variogram  results 
used  to  develop  these  figures  are  presented  in  Appendix  D.  The  variograms  of 
Figures  15  through  20  are  representative  of  the  types  of  results  observed  for 
variograms  in  this  study. 

25.  Many  of  the  cross  sections  of  the  DS  +  122  joint  surface  were  asso¬ 
ciated  with  variograms  showing  a  behavior  attributable  to  the  data  being  non¬ 
stationary  (see  Appendix  B  for  an  explanation  of  the  Intrinsic  Hypothesis) . 
This  type  of  beha%Tior  results  when  a  change  in  data  values  occurs  in  a  regular 
fashion  over  a  particular  direction,  such  as  elevation  values  which  increase 
as  one  travels  upslope.  Cross  sections  which  yielded  such  variograms  did  show 
a  distinct  inclination.  Non-stationary  variograms  are  parabolic  and  concave 
upward  in  shape.  Figure  15  shows  a  variogram  of  this  type  for  cross  section 
BG.  In  Figure  13,  this  cross  section  is  seen  to  be  inclined. 

26.  For  this  study,  and  non-stationary  variogram  in  general,  one  proce¬ 
dure  to  make  the  data  more  amenable  to  geostatistical  analysis  is  to  filter 
the  non-stationary  component  from  the  data.  If  the  data  are  denoted  as  Z(x), 
two  components  comprise  these  data,  S(x)  and  N(x),  such  that 

Z(x)  =  S(x)  +  N(x)  (5) 

where  S(x)  is  the  stationary  component  and  N(x)  is  the  non-stationary  compo¬ 
nent.  For  this  study,  S(x)  is  the  elevation  value  due  to  surface  asperities 
and  N(x)  is  the  elevation  value  due  to  the  inclination  of  the  cross  section. 

27.  For  an  enhanced  analysis,  N(x)  was  filtered  from  the  cross  section 
data,  such  as  cross  section  BG,  to  yield 

S (x)  =  Z(x)  -  N(x)  (5a) 

The  procedure  used  to  remove  N(x)  is  shown  in  Figure  21.  The  angle  of  incli¬ 
nation  was  calculated,  and  the  elevation  of  the  plane  defined  by  this  angle 
was  removed  from  each  of  the  cross  section  elevation  values.  This  resulted  in 
filtered  cross  sections.  A  variogram  computed  for  cross  section  BG,  after 
filtering,  is  shown  in  Figure  22. 
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Figure  15.  Variogram  for  cross  section  BG,  an  example  of 
non-stationary  behavior  in  the  variogram 
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Figure  16.  Variogram  for  cross  section  AF.  Variogram  is  stationary, 
but  at  h  =  400  cm  ,  it  becomes  non-stationary 
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Figure  17.  Variogram  for  cross  section  AL,  an  example  of 

stationary  variogram 
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Figure  18.  Variogram  for  cross  section  AS',  an  example  of 

stationary  behavior 
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Figure  19.  Variogram  for  cross  section  BA,  an  example  of 

stationary  behavior 
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Figure  20.  Variogram  for  cross  section  BB,  an  example  of  a 
non-stationary  variogram 
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Figure  21.  Diagram  of  the  procedure  used  Lo  remove 
the  influence  of  the  plane  of  inclination  from 
cross-sectional  data 
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Figure  22.  Variogram  for  cross  section  BG  subsequent  to 
the  removal  of  the  inclination 
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Results  of  the  Variogram  Analysis 


28.  Stationary  variograms  (Figures  17,  18,  19,  and  22)  show  that  asper¬ 
ity  undulations  along  a  joint  surface  are  spatially  correlated.  Table  4  lists 
the  properties  of  these  four  variograms.  The  sill,  a  more  or  less  constant 
value  the  variogram  attains  at  a  certain  separation  distance,  of  each  vario¬ 
gram  is  slightly  larger  than  the  variance  of  the  elevation  values  along  each 
cross  section.  The  range,  that  separation  distance  at  which  the  variogram 
attains  the  sill,  of  each  variogram  shows  the  spatial  continuity  of  asperity 
values;  there  is  not  spatial  correlation  for  separation  distances  greater  than 
the  range.  The  nugget  value,  the  intercept  of  the  variogram  with  the  ordinant 
axis,  shows  the  level  of  noise  (random  variations)  present  in  the  elevation 
data;  for  each  cross  section,  the  nugget  value  is  small  as  expected  for  this 
type  of  precisely  measured  data.  Finally,  each  variogram  is  compared  to  the 
fractal  dimension  for  a  particular  cross  section.  There  appears,  on  the  basis 
of  at  least  these  four  variograms,  to  be  little  correlation  between  the  vario¬ 
grams  and  the  fractal  dimension.  The  fractal  dimension  offers  a  unique  method 
for  the  description  of  surface  roughness.  The  variogram,  while  yielding  no 
clear  information  on  surface  roughness,  analyzes  other  features  of  a  joint 
surface,  particularly  the  spatial  distribution  of  asperities. 

Table  4 

Properties  of  Stationary  Variograms 


Cross  Section 

Nugget 

Sill 

Range  (cm) 

Fractal 

Dimension 

AL 

0.0 

1400.00 

300.00 

1.0100 

AS' 

0.0 

650.00 

330.00 

1.0052 

BA 

0.0 

130.00 

150.00 

1.0000 

BG  (filtered) 

0.0 

6.00 

80.00 

1.0000 

Potential  Application 

of  the 

Variogram:  Kriging 

29.  One  value  result  of  this  study  was  the  demonstration  that  vario¬ 
grams  can  be  developed  for  elevation  data  measured  along  a  joint  surface. 
Because  these  data  are  shown  to  be  spatially  correlated,  this  correlation  can 
be  used  to  estimate  elevation  values  at  locations  where  such  values  were  not 
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measured.  This  estimation  process  is  achieved  using  the  geostatistical  esti¬ 
mator  known  as  kriging.  This  technique  is  described  in  Appendix  B. 

30.  For  geostatistical  analyses,  such  as  were  reviewed  in  this  study, 
kriging  is  useful  in  mapping  available  data.  For  the  cross  sections  used  for 
this  study,  kriging  could  be  used  to  develop  greater  detail  in  these  plots. 

For  example,  kriging  could  be  used  to  estimate  an  elevation  value  between  two 
measured  locations.  This  estimate  would  be  based  on  the  two  known  elevation 
measurements  and  the  variogram.  In  this  manner,  the  spatial  correlation  and 
continuity  shown  by  the  asperities  could  be  used  to  better  define  each 
asperity  (i.e.,  to  fill  in  data  between  measurement  points). 

31.  More  importantly,  this  study  relied  on  a  plotted  contour  map  of  the 
DS  +  122  joint  surface  as  a  source  of  data.  Elevation  values  of  such  detail, 
already  mapped,  are  rare.  More  often,  elevation  values  are  typically  measured 
across  a  joint  surface  at  irregular  locations  with  a  rather  sparse  sampling 
interval.  The  results  of  this  study  imply  that  a  variogram  could  be  developed 
for  this  type  of  sparse  data.  Kriging  could  be  used  to  estimate  elevation 
values  at  intersections  of  a  regular  grid  superimposed  over  the  joint  surface. 
In  this  manner,  kriging  could  be  useful  in  developing  the  type  of  contour  map 
which  served  as  the  data  base  for  this  study  without  the  need  for  photogramme- 
trically  mapping  the  entire  surface. 

32.  Kriging  is  often  used  for  gridding  purposes  (i.e.,  the  interpola¬ 
tion  of  data  at  intersections  of  a  regular  grid;  Olea  1974).  This  is  done 
especially  in  preparation  for  data  contouring;  an  algorithm  for  contouring 
usually  requires  a  regular  grid  of  data.  A  rock  joint  is  a  two  dimensional 
surface  and  can  be  represented  by  a  two  dimensional  grid.  This  report,  how¬ 
ever,  examines  one  dimensional  profiles  along  certain  directions  across  a 
joint  surface,  Kriging  might  be  used  in  a  broader  investigation  to  form  a  two 
dimensional  digital  grid  of  a  joint  surface.  Such  a  grid  could  then  be  used 
to  calculate  a  fractal  dimension  for  the  entire  surface,  rather  than  for 
discrete  profiles.  This  might  lead  to  a  significant  contribution  to  the  pro¬ 
cedures  used  for  the  analysis  of  sliding  stability  of  rock  masses. 
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PART  IV:  COMPARISON  OF  THE  USE  OF  THE  FRACTAL  DIMENSION 
AND  THE  VARIOGRAM  TO  ACCEPTED  ROCK  SURFACE 
ROUGHNESS  EVALUATION  TECHNIQUES 

33.  It  is  emphasized  that  the  use  of  the  fractal  dimension  and  the 
variogram  supplements,  but  does  not  replace,  classic  surface  roughness  evalu¬ 
ation  techniques.  It  has  been  noted  in  the  past  that  different  types  of  joint 
and  fault  surfaces  display  quite  different  surface  roughness  patterns  (Patton 
1966).  It  was  further  observed  that  the  surface  roughness  varies  with  direc¬ 
tion  across  a  rock  surface  (Patton  and  Deere  1970).  These  orientation  differ¬ 
ences  are  especially  important  for  sliding  stability  analyses. 

34.  Surface  roughness  is  caused  by  irregularities,  also  called  asperi¬ 
ties.  There  are  usually  several  orders  of  magnitudes  of  surface  asperities, 
from  microscopic  pits  to  megascopic  folds  (Piteau  1970).  For  the  analysis  of 
sliding  stability,  it  is  often  assumed  that  small  magnitude  asperities  bear 
the  initial  shear  load,  but  quickly  shear.  Subsequent  to  this,  the  megascopic 
asperities  bear  the  load  and  control  the  sliding  stability  (Hoek  and  Bray 
1974). 

35.  The  different  magnitudes  of  asperities  are  easily  recognized  in 
Figure  23,  which  is  reproduced  from  Patton  and  Deere  (1970).  This  figure 
shows  an  actual  cross  section  along  a  bedding  plane  and  is  analogous  to  the 


Figure  23.  Cross  section  form  Patton 
and  Deere  (1970) 
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cross  sections  plotted  for  the  DS  +  122  joint  surface.  The  fractal  dimension 
of  the  cross  section  of  Figure  23  is  1.01. 

36.  It  is  evident  from  this  example  that  the  fractal  dimension  does  not 
clearly  identify  the  different  orders  of  asperity  magnitudes.  The  fact  that 
there  are  different  magnitudes  of  asperities,  however,  does  follow  the  concept 
of  self-similarity  (Mandelbrot  1982).  Each  of  the  small  asperities  is  self¬ 
similar  to  the  next  larger  size  of  asperity,  and  so  on.  It  is  hypothesized 
that  the  more  magnitudes  of  asperity  sizes  a  joint  surface  has,  the  larger 
will  be  its  fractal  dimension.  In  this  sense,  the  fractal  dimension  might 
serve  as  an  indicator  of  the  presence  of  different  magnitudes  of  asperities. 

37.  The  variogram  is  more  valuable  for  the  analysis  of  asperity  magni¬ 
tudes.  The  cross  section  of  Figure  23  is  replotted  in  Figure  24  to  remove  the 
inclination  and  increase  the  scale.  Assigned  for  this  study  to  be  consistent 
with  the  horizontal  scale,  an  elevation  scale  is  shown.  This  scale  was  not 
originally  provided  by  Patton  and  Deere  (1970),  but  can  be  assumed  to  be  cor¬ 
rect  based  on  the  inclination  values  noted  in  Figure  23.  A  variogram  was 
developed  for  this  cross  section  and  is  plotted  in  Figure  25. 

38.  Of  particular  interest  to  this  study,  the  variogram  of  Figure  25 
shows  a  periodic  behavior.  The  variogram  increases  initially  away  from  the 
origin,  attains  a  plateau  (the  sill),  then  decreases,  but  subsequently 
increases  to  a  second  plateau.  The  distance,  h  ,  between  the  origin  and  the 
point  at  which  the  variogram  begins  Its  second  increase  is  approximately 

100  cm. 

39.  In  Figure  24,  there  is  a  100  cm  distance  between  the  peaks  of  the 
megascopic  asperities.  The  variogram  of  Figure  25  is  interpreted  to  represent 
the  spatial  correlation  between  the  small  asperities  (the  first  increase  and 
plateau)  and  the  spatial  correlation  between  the  larger  and  more  widely  spaced 
asperities  (the  second  increase  and  plateau).  This  is  a  valuable  application 
of  the  variogram  because  it  appears  from  this  example  that  the  variogram  can 
help  to  describe  the  various  magnitudes  of  asperities. 

40.  For  the  cross  sections  of  the  DS  +  122  joint  surface,  several 
variograms  displayed  the  type  of  behavior  shown  by  the  variogram  of  Figure  25. 
Two  of  these  variograms  are  again  plotted  in  Figures  26  and  27;  the  cross  sec¬ 
tions  are  plotted  in  Figure  28.  For  cross  section  AT.,  the  variogram  identi¬ 
fies  a  second,  larger  asperity  group  with  an  average  spacing  of  600  cm.  For 
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Figure  26.  Variogram  for  cross  section  AL,  presented  here  to 
emphasize  its  periodic  behavior 
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Figure  27.  Variogram  for  cross  section  AS',  presented  here 
to  emphasize  its  subtle  periodic  behavior 
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Figure  28.  Cross  sections  AL  and  AS',  plotted  for  comparison 
with  the  variograms  of  Figures  26  and  27 

cross  section  AS',  the  spacing  of  the  larger  asperities  is  found  from  the 
variogram  to  be  680  cm. 

41.  Directional  characteristics  of  surface  roughness,  as  noted  by 
Patton  and  Deere  (1970),  can  be  evaluated  using  the  fractal  dimension  or  the 
variogram.  Cross  sections  along  the  DS  +  122  joint  surface  had  specific  ori¬ 
entations.  Hence,  each  computed  fractal  dimension  and  variogram  indicated  a 
specific  directional  characteristic.  The  combination  of  fractal  dimension  and 
variogram  is  ideally  suited  for  the  description  of  joint  properties.  The 
qualitative  descriptions  offered  by  Patton  (1966),  Patton  and  Deere  (1970), 
and  Piteau  (1970)  can  be  quantitatively  described  by  the  fractal  dimension  and 
the  variogram. 
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PART  V:  APPLICATION  OF  THE  FRACTAL  DIMENSION  AND  VARIOGRAM 
TO  STRING  LINE  DATA  FOR  ROCK  SURFACES  NEAR 
LIBBY  DAM,  MT 

42.  Previous  sections  of  this  report  reviewed  the  application  of  the 
fractal  dimension  and  variogram  to  the  analysis  of  the  DS  +  122  slide  surface 
adjacent  to  Libby  Dam,  MT.  Subsequent  to  this  slide,  ground-based  reconnais¬ 
sance  was  made  of  older  joint  and  bedding  plane  surfaces  in  proximity  to  the 
DS  +  122  slide.  During  the  reconnaissance,  string  line  measurements  were  made 
of  several  rock  surfaces.  This  process  involved  the  placement  of  a  string 
approximately  parallel  to  and  above  the  rock  surface.  Elevations  were  mea¬ 
sured  from  the  string  line  to  the  rock  surface  at  0.5  ft*  intervals.  The  rock 
surfaces  for  which  string  line  data  were  collected  are  listed  in  Table  5. 

These  data  are  plotted  in  Figures  29  through  33. 

43.  These  string  line  data  provided  a  different  means  for  evaluating 
the  fractal  dimension  and  spatial  continuity  of  rock  surfaces.  The  data  used 
previously  to  describe  the  fractal  dimension  of  the  DS  +  122  slide  were 
derived  from  a  highly  detailed  contour  map  of  the  rock  surface  plotted  at  2  cm 
contour  intervals.  At  times,  this  contour  map  was  difficult  to  read.  More¬ 
over,  the  map  was  plotted  from  stereo  photo  pairs,  and  the  true  data  resolu¬ 
tion  is  estimated  by  observations  during  this  study  to  be  20  cm.  Hence,  the 
string  line  data  may  be  a  higher  resolution  data  source  and  possibly  more 
reliable  than  the  data  obtained  from  the  contour  map. 

Fractal  Dimension  of  String  Line  Cross  Sections 

44.  Fractal  dimensions  were  calculated  for  the  rock  surfaces  listed  in 
Table  5.  The  procedure  used  was  the  one  reviewed  for  the  DS  +  122  cross  sec¬ 
tion  using  the  technique  whereby  the  fractional  remainder  is  added  to  the 
number  of  counts,  N  ,  for  each  segment  length,  y  .  The  results  are  listed 
in  Table  6. 

45.  For  this  application,  log^CN)  versus  log^(y)  plots  were  not  drawn 
to  determine  the  fractal  dimension.  Instead,  a  least  squares  regression 


*  A  table  of  factors  for  converting  non-SI  units  of  measurement  to  SI 
(metric)  units  is  presented  on  page  4. 
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Table  5 

List  of  Rock  Surfaces  Associated  with  String  Line  Data 


Rock 

Surface  (Names  from  Nearest 

Topographic  Feature) 

Typt  ui  Feature 

1. 

Backus  Notch  No.  1 

Bedding  plane 

2. 

Backus  Creek  Notch  No.  1 

Bedding  plane 

3. 

DS  +  122,  C 

Joint. 

4. 

Island  Notch 

Joint 

5 . 

DS  +  122  (minor) 

Joint 

6. 

DS  +  122,  No.  2 

Joint 

7. 

914  Rib* 

Joint 

8. 

Dunn  Creek  Notch 

Joint 

9. 

930  joint 

Joint 

10. 

Island  Notch 

Joint/bedding  plane 
intersection 

11. 

Island  Notch 

Bedding  plane 

12. 

Wolf  Creek  Jet. 

Joint 

13. 

Dunn  Creek  Notch 

Bedding  plane 

14. 

Old  Notch  No.  1 

Joint 

15. 

Backus  Notch  No.  1 

Joint 

*  Note:  The  914  Rib  joint,  the  930  joint  and  the  Old  Notch  No.  1  joint  are 
prehistoric  slide  scars.  The  DS  +  122  C,  DS  +  122  minor  and  DS  +  122  No.  2 
profiles  are  in  situ  traces  of  joints  exposed  on  the  primary  DS  +  122 
bedding  surface  and  are  members  of  the  same  set  as  the  failed  joint.  The 
remainder  of  the  surfaces  are  natural  surfaces  exposed  during  construction 
of  the  highway  and  the  Burlington  Northern  railroad  relocations. 
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Figure  29.  Plot  of  string  line  cross  sections:  Backus 
Notch  No.  1  (bedding);  Backus  Creek  Notch  No.  1  (bedding); 
and  DS  +  122  (C  joint) 
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Figure  30.  Plot  of  string  line  cross  sections: 
Island  Notch  (joint);  DS  +  122  (minor  joint); 
and  DS  +  122  (joint  No.  2) 
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Figure  31.  Plot  of  string  line  cross  sections:  914  Rib  (joint); 
Dunn  Creek  Notch  (joint);  and  930  (joint) 
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Figure  32.  Plot  of  string  line  cross-sections:  Island 
Notch  (joint/bedding  plane  intersection);  Island  Notch 
(bedding);  and  Wolf  Creek  Jet.  (joint) 
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Figure  33.  Plot  of  string  line  cross  sections:  Dunn  Creek 
Notch  (bedding);  Old  Notch  No.  1  (joint);  and  Backus  Notch 

No.  1  (joint) 
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Table  6 

Fractal  Dimension  Calculation  for  String  Line  Cross  Section 


Number  of 

Counts,  N 

Fractal 

Rock  Surface 

*y  =  2.0 

y  -i.5 

y  =  1.0 

y  =  0.5 

Dimension 

Backus  Notch  No.  1  (bed) 

46.375 

61.847 

92.790 

185.84 

1.001320 

Backus  Creek  Notch  (bed) 

31.835 

42.493 

63.615 

127.80 

1.002256 

DS  +  122  (C  joint) 

50.330 

67.093 

100.73 

201.99 

1.002557 

Island  Notch  (joint) 

10.265 

13.690 

20.580 

40.210 

1.002831 

DS  +  122  (minor  joint) 

50.070 

66.860 

100.38 

201.13 

1.002907 

DS  +  122  (joint  No.  2) 

50.175 

66.893 

100.52 

201.52 

1.003175 

914  Rib  (joint) 

33.080 

44.173 

66.370 

132.92 

1.003190 

Dunn  Creek  Notch  (joint) 

10.020 

13.400 

20.050 

40.300 

1.003220 

930  (joint) 

35.025 

46.720 

70.280 

141.08 

1.005294 

Island  Notch  (joint/bed) 

9.290 

12.380 

18.680 

37.420 

1.005762 

Island  Notch  (bed) 

42.625 

56.840 

85.280 

171.970 

1.006291 

Wolf  Creek  Jet.  (joint) 

37.750 

50.490 

76.180 

152.84 

1.008842 

Dunn  Creek  Notch  (bed) 

14.165 

18.990 

28.400 

57.560 

1.010249 

Old  Notch  No.  1  (joint) 

38.305 

51.487 

77.700 

157.20 

1  .017850 

Backus  Notch  No.  1  (joint) 

41.230 

54.980 

82.730 

166.08 

1.005408 

■'or  y  .re  in  feet 

;  these 

correspond 

to  60.96  cm, 

45.72  cm 

9 

and  IS. 24  cm. 


<-v  ,  r  ;  vis  used  to  calculate  the  fractal  dimension.  A  linear  regression 

ccua-nn  i  a  ess  oined  of  the  form: 


log10(N)  =  D-q  +  Ojlog^Cy)  (6) 

In  this  equation,  is  the  value  of  log^(N)  when  log^(y)  is  zero  and 

is  the  slope  of  this  line.  Because  the  fractal  dimension,  D  ,  is  the 
negative  of  the  slope  of  this  line,  then 

D  =  -D1  (7) 
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46.  Using  the  principle  of  least  squares  regression,  D  is  calculated 
as  (McCuen  1985,  p  187): 

E(Log10(N)Log10(y))  -  (ZLog1()(N)ZLog10(y))/K 

D  =  -D.  =  -  - - - x -  (8) 

Z(Log10(y)r  -  (ZLog10(y))  /K 

where  K  is  the  number  of  different  trial  segment  lengths,  y  ,  used  to  cal¬ 
culate  D  .  In  Table  6,  for  example  K  is  equal  to  4  and  the  fractal  dimen¬ 
sion  values  listed  were  calculated  using  equation  (8).  The  fractal  dimension 
values  range  from  1.000132  for  the  Backus  Notch  No.  1  bedding  surface  to 
1.01785  for  the  Old  Notch  No.  1  joint  surface.  These  values  are  comparable  to 
those  determined  for  the  DS  +  122  joint  surface  listed  in  Table  2  for  the 
Method  2  calculation. 


Computer  Calculation  of  the  Fractal  Dimension 


47.  Results  presented  in  Tables  2  and  6  were  developed  by  selecting  a 
segment  length,  y  ,  then  physically  walking  a  pair  of  dividers  opened  a  dis¬ 
tance  equal  to  y  along  a  plotted  cross  section  or  string  line.  This  is  a 
tedious  task  and  is  susceptible  to  error  during  the  counting  of  the  number  of 
steps,  N  ;  the  calculation  of  the  remainder  at  the  end  of  the  plot;  or  from 
stretching  of  the  paper  on  which  the  cross  sections  are  plotted  due  to  humid¬ 
ity;  etc.  A  computer  program  was  developed  to  enable  easier  calculation  of 
the  fractal  dimension  by  a  numerical  algorithm  of  counting  segment  lengths  in 
addition  to  equation  (8).  This  program  is  presented  in  Appendix  E. 

48.  Input  to  the  program  consists  of  the  number  of  segment  lengths, 

y  ,  to  be  considered;  specification  of  these  segment  lengths;  and  the  eleva¬ 
tion  and  distance  coordinates  along  a  cross  section  or  string  line.  This  pro¬ 
gram  is  limited  at  present  to  one  dimensional  analyses.  The  algorithm  on 
which  the  program  is  based  is  graphically  presented  in  Figure  34. 

49.  Fractal  dimensions  calculated  for  the  rock  surfaces  listed  in 
Table  5  using  the  computer  program  of  Appendix  E  are  listed  in  Table  7. 
Although  the  fractal  dimension  values  do  not  appear  to  be  much  different  than 
those  shown  in  Table  6,  a  different  result  was  obtained  when  analyzing  the 
cross  sections  in  terms  of  roughness.  The  fractal  dimension  values  are  so 
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Table  7 

Computer  Calculation  of  String  Line  Fractal  Dimensions 


Number  of 

Counts,  N 

Fractal 

Rock  Surface 

*y  =  2.0 

y  =  1.5 

y  =  1.0 

y  =  0.5 

Dimension 

Backus  Notch  No.  1  (bed) 

46.303 

61.755 

92.682 

185.48 

1.001059 

Backus  Creek  Notch  (bed) 

31.785 

42.391 

63.603 

127.30 

1.000891 

DS  +  122  (C  joint) 

50.334 

67.105 

100.78 

202.14 

1.003045 

Island  Notch  (joint) 

10.252 

13.670 

20.504 

41.013 

1.000107 

DS  +  122  (minor  joint) 

50.089 

66.824 

100.337 

201.28 

1.003369 

DS  +  122  (joint  No.  2) 

50.127 

66.887 

100.55 

201.44 

1.003469 

914  Rib  (joint) 

33.048 

44.171 

66.248 

132.88 

1.003334 

Dunn  Creek  Notch  (joint) 

10.028 

13.387 

20.110 

40.377 

1.004735 

930  (joint) 

35.018 

46.731 

70.256 

140.93 

1.004516 

Island  Notch  (joint/bed) 

9.278 

12.389 

18.586 

37.285 

1.003104 

Island  Notch  (bed) 

42.592 

56.829 

85.272 

170.72 

1.001409 

Wolf  Creek  Jet.  (joint) 

37.850 

50.555 

76.094 

152.57 

1.005641 

Dunn  Creek  Notch  (bed) 

14.189 

19.076 

28.458 

57.588 

1.008608 

Old  Notch  No.  1  (joint) 

38.345 

51.624 

77.055 

158.28 

1.021994 

Backus  Notch  No.  1  (joint) 

41.231 

55.004 

82.666 

166.86 

1.004233 

*  Units  for  y  are  in  feet 

;  these 

correspond 

to  60.96  cm, 

45.72  cm 

» 

30.48  cm,  and  15.24  cm. 


close  in  value  to  1  that  considerable  error  can  occur  when  measuring  these 
values  by  hand.  Small  changes  in  the  number  of  counts,  N  ,  determined  for  a 
segment  length,  y  ,  can  result  in  a  significantly  different  fractal  dimen¬ 
sion.  A  computer  calculation  of  the  fractal  dimension  is,  therefore,  pre¬ 
ferred  over  a  hand  calculation  procedure. 

Variogram  Analysis  of  String  Line  Cross  Sections 

50.  Variograms  were  calculated  for  each  rock  surface  listed  in  Table  5. 
These  variograms  are  presented  in  Appendix  F.  Table  8  lists  the  string  line 
cross  sections,  their  fractal  dimensions,  and  variogram  characteristics.  The 
following  observations  are  made  of  these  variograms: 

a.  The  variogram  nugget  is  strongly  correlated  with  the  fractal 
dimension  (the  correlation  coefficient  is  0.979). 

b.  The  variogram  sill  is  strongly  correlated  with  the  fractal 
dimension  (the  correlation  coefficient  is  0.914). 
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Table  8 

Varlogram  Results  for  String  Line  Cross  Sections 


Rock  Surface 

Fractal  . 

sjc 

Dimension 

Nugget 

Variogram 
Sill  Range 

Behavior 

Backus  Notch  No.  1  (bed) 

1.001059 

0.0009 

0.0149 

10.00 

2 

periods 

Backus  Creek  Notch  (bed) 

1.000891 

0.0008 

0.0848 

27.00 

0 

periods 

DS  +  122  (C  joint) 

1.003045 

0.0039 

0.0450 

16.00 

4 

periods 

Island  Notch  (joint) 

1.000107 

0.0001 

0.0010 

7.20 

1 

period 

DS  +  122  (minor  joint) 

1.003369 

0.0026 

0.0651 

16.00 

0 

periods 

DS  +  122  (joint  No.  2) 

1.003469 

0.0027 

0.1590 

25.00 

0 

periods 

914  Rib  (joint) 

1.003334 

0.0028 

0.0200 

30.00 

2 

periods 

Dunn  Creek  Notch  (joint) 

1.004735 

0.0035 

0.0200 

6.50 

2 

periods 

930  (joint) 

1.004516 

0.0050 

0.2680 

14.00 

0 

periods 

Island  Notch  (joint/bed) 

1.003104 

0.0025 

0.0552 

6.00 

0 

periods 

Island  Notch  (bed) 

1.001409 

0.0016 

0.0164 

23.00 

2 

periods 

Wolf  Creek  Jet.  (joint) 

1  .005641 

0.0040 

0.0349 

14.00 

0 

periods 

Dunn  Creek  Notch  (bed) 

1.008608 

0.0134 

0.0677 

5.50 

2 

periods 

Old  Notch  No.  1  (joint) 

1.021994 

0.0254 

1.8800 

12.20 

0 

periods 

Backus  Notch  No.  1  (joint) 

1.004233 

0.0045 

0.3840 

16.20 

0 

periods 

Correlation  Coefficient 


1.  Fractal  Dimension  and  Nugget  0.979 

2.  Fractal  Dimension  and  Sill  0.914 

3.  Fractal  Dimension  and  Range  -0.217 

4.  Fractal  Dimension  and  Behavior  -0.206 


*  The  fractal  dimension  values  shown  are  from  Table  1 .  Units  for  the  range 
are  feet.  The  behavior  category  describes  the  appearance  of  the  variogram 
relative  to  Figure  25. 


c.  There  is  a  poor  correlation  between  variogram  periodicity  (a 
behavior  similar  to  Figure  25)  and  the  fractal  dimension;  the 
correlation  coefficient  is  -0.206;  hence,  the  variogram  is  less 
likely  to  be  periodic  for  large  fractal  dimensions  (a  negative 
correlation) . 

d.  The  variogram  range  is  poorly  correlated  with  the  fractal 
dimension  (the  correlation  coefficient  is  -0.217). 

From  these  variogram  analyses,  the  strongest  correlation  is  found  between  the 

variogram  nugget  (an  indicator  of  data  randomness)  and  the  fractal  dimension. 

A  strong  correlation  is  also  found  between  the  variogram  sill  and  the  fractal 

dimension.  Other  variogram  characteristics  are  seemingly  unrelated  to  the 

fractal  dimension. 

51.  A  numerical  example  is  forwarded  to  explain  the  relationship 
between  the  fractal  dimension  and  variogram  nugget.  A  random  number  generator 
was  used  to  simulate  elevation  values  at  a  spacing  of  0.5  ft  along  a  ficti¬ 
tious  50  ft  string  line.  These  simulated  elevation  values  were  scaled  to  have 
a  range  between  0  and  4  ft.  The  result  is  plotted  in  Figure  35.  Using  the 
computer  program  listed  in  Appendix  E,  the  fractal  dimension  was  calculated  to 
be  1.2952.  The  variogram  of  this  cross  section,  shown  in  Figure  36,  is  ran¬ 
dom,  as  expected  for  random  numbers;  the  nugget  value  is  equal  to  sample 
variance,  1.20. 

52.  Three  different  low  pass  filters  were  applied  to  these  random  ele¬ 
vation  values.  These  filtered  results  are  plotted  in  Figure  37,  and  vario- 
grams  are  plotted  in  Figures  38  through  40.  The  filtering  process  results  in 
spatially  correlated  elevation  values;  hence,  the  variograms  develop  an 
appearance  similar  to  the  spherical  shape  (Appendix  B) .  As  the  filter  size 
increases,  the  variance  of  the  data  decreases;  thus,  the  variogram  sill 
decreases.  This  smoothing  results  in  a  lower  fractal  dimension  and  a  reduc¬ 
tion  in  data  randomness;  hence,  the  nugget  value  decreases.  This  demonstrates 
that  the  nugget  and  sill  values  of  a  variogram  are  directly  proportional  to 
the  fractal  dimension;  as  the  fractal  dimension  approaches  unity,  the  vario¬ 
gram  nugget  and  sill  approach  zero. 
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Figvire  35.  Artificial  string  line  cross  section  formed  using  random 

elevation  values 


gamma  (h) 
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•  111 

8  12  16  20 

distance,  h 

Figure  36. 

Variogram 

for  string  line  cross  section  shown  in 

Figure  35 
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Figure  37.  Three  filtered  results  developed  from  the  string  line 
cross  section  shown  in  Figure  35;  filter  sizes  are  shown 
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Figure  38.  Variogram  for  the  3X1  filtered  string  line  cross 

section 
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Figure  39.  Variogram  for  the  7X1  filtered  string  line  cross 

section 
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Figure  40.  Variogram  for  the  15X1  filtered  string  line  cross 

section 
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PART  VI:  REGRESSION:  THE  RELATIONSHIP  BETWEEN  THE  FRACTAL 
DIMENSION  AND  JOINT  ROUGHNESS  COEFFICIENT 


53.  Joint  roughness  (and  the  roughness  of  bedding  planes)  is  a  factor 
in  the  sliding  stability  of  rock  masses.  The  rougher  the  surface  of  a  joint 
is,  the  greater  stability  it  has  against  sliding.  An  accepted  value  which 
describes  the  roughness  of  rock  surfaces  is  known  as  the  Joint  Roughness 
Coefficient  or  JRC  (Barton  and  Choubey  1977). 

54.  Typical  roughness  profiles  for  various  JRC  values  are  shown  in  Fig¬ 
ure  41,  reproduced  from  Barton  and  Choubey  (1977).  JRC  values  were  assigned 
to  each  cross  section  listed  in  Table  5  using  Figure  41.  This  assignment  is 


Figure  41.  JRC  profiles  (Barton  and  Choubey 

1977) 
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shown  in  Table  9.  It  is  noted  that  one  cross  section  is  assigned  a  value  of 
22,  which  is  larger  than  any  profile  shown  by  Barton  and  Choubey  (1977).  This 
high  value  was  extrapolated  based  on  the  nature  of  the  other  Libby  Dam  area 
sections.  Moreover,  the  sections  shown  by  Barton  and  Choubey  (1977)  refer  to 
the  laboratory  specimens  they  examined  and  are  not  thought  to  represent  maxi¬ 
mum  possible  values.  Each  cross  section  is  listed  along  with  the  assigned  JRC 
value,  the  fractal  dimension  and  the  variogram  nugget.  A  fictitious,  ideal 
cross  section  is  also  included  for  the  subsequent  discussion.  This  fictitious 
rock  surface  is  an  ideal  plane  with  a  JRC  value  equal  to  zero,  a  fractal 
dimension  equal  to  one,  and  a  nugget  equal  to  zero.  The  purpose  of  including 
this  ideal  surface  is  explained  subsequently. 

55.  From  the  information  listed  in  Table  9,  a  regression  equation  can 
be  developed  to  predict  a  JRC  value  based  on: 

a.  The  fractal  dimension  alone. 

b.  The  variogram  nugget  alone. 

c.  The  fractal  dimension  and  variogram  nugget. 

The  development  of  such  a  regression  equation  has  a  simple  objective.  It  is  a 
difficult  exercise  in  judgment  to  use  Figure  41  to  assign  a  JRC  value  to  a 
specific  cross  section.  By  developing  a  regression  of  JRC  on  the  fractal 
dimension  and/or  the  variogram  nugget,  the  JRC  value  can  be  estimated  rather 
than  visually  assigned. 

56.  The  calculation  of  JRC  is  not  a  new  approach.  A  calculation  proce¬ 
dure  was  developed  based  on  laboratory  shear  tests  of  rock  samples  (Barton  and 
Choubey  1977) : 


where 


arctant  (i/cr^)  - 
LogI0(JCS)/„N) 


(9) 


t  =  shear  stress 

o..  =  effective  normal  stress 
N 

<J>k  =  basic  friction  angle 
JCS  =  joint  wall  compression  strength 

The  solution  to  equation  (9)  requires  the  use  of  laboratory  rock  shear  tests. 

57.  The  development  of  a  regression  equation  based  on  the  fractal 
dimension  and  variogram  nugget  requires  only  elevation  data  measured  along  a 


Table  9 

Joint  Roughness  Coefficient  (JRC)  for  String  Line  Cros3  Section 


Rock  Surface 

JRC* 

D 

CO 

Fictitious  Ideal 

0 

1.0000 

0.0000 

Backus  Notch  No.  1  (bedding) 

3 

1.0011 

0.0009 

Backus  Creek  Notch  No.  1  (bedding) 

2 

1.0009 

0.0008 

DS  +  122  (C  joint) 

4 

1.0031 

0.0039 

Island  Notch  (joint) 

1 

1.0001 

0.0001 

DS  +  122  (minor  joint) 

4 

1.0034 

0.0026 

DS  +  122  (joint  No.  2) 

5 

1.0035 

0.0027 

914  Rib  (joint) 

3 

1.0033 

0.0028 

Dunn  Creek  Notch  (joint) 

7 

1.0047 

0.0035 

930  (joint) 

8 

1.0045 

0.0050 

Island  Notch  (joint/bed) 

5 

1.0031 

0.0025 

Island  Notch  (bedding) 

2 

1.0014 

0.0016 

Wolf  Creek  Jet.  (joint) 

8 

1.0056 

0.0040 

Dunn  Creek  Notch  (bedding) 

13 

1.0086 

0.0134 

Old  Notch  No.  1  (joint) 

22 

1.0220 

0.0254 

Backus  Notch  No.  1  (joint) 

6 

1.0042 

0.0045 

*  JRC  is  the  Joint  Roughness 

Coefficient;  D 

is  the 

fractal  dimension; 

and  CO 

is  the  variogram  nugget.  For  this  table, 

N  is 

16. 

consistent  linear  direction. 

These  data  can 

be  collected  in  the 

field 

without 

the  need  for  expensive  laboratory  equipment. 

The  development  of 

such 

an  equa- 

t^on  will  not  duplicate,  but 

will  supplement 

or  be 

an  alternative 

to 

equation  (9) . 

Regression  Analysis 

58.  Regression  analysis  involves  the 

estimation  of  a  dependent 

variable 

using  one  or  several  independent  variables. 

Such  a 

procedure  begins  by  plott- 

ing  the  dependent  variable  versus  the  independent  variables,  one  plot  per 
independent  variable.  Such  a  plot  reveals  the  type  of  correlation  between  the 


dependent  and  independent  variables.  There  may  be  a  lack  of  correlation,  or 
the  relationship  may  be  linear,  quadratic  or  some  other  functional  form. 

59.  Suppose,  for  example,  the  two  variables  are  linearly  correlated  and 
only  one  independent  variable  is  considered.  Figure  42  is  an  example.  In 
this  case,  the  dependent  variable,  y  ,  can  be  estimated  from  the  independent 
variable,  x  ,  as 

y  =  AQ  +  Alx  (10) 

Once  the  form  of  this  equation  has  been  chosen,  it  remains  to  solve  for  the 
coefficients,  A^  and  A^  ,  to  define  the  equation  explicitly. 

60.  In  regression  analysis,  one  possible  method  for  the  solution  of  the 
coefficients  is  the  method  of  least  squares.  This  involves  the  minimization 


Variable  x 


Figure  42.  Example  case  of  two  linearly  correlated 

variables 


of  the  errors  when  estimating  the  dependent  variable,  y  .  To  explain  this 
procedure,  Figure  42  is  recalled.  Each  point  plotted  on  this  figure  is  asso¬ 
ciated  with  two  known  values,  x  and  y  .  At  each  of  these  points,  an  esti¬ 
mate  of  y  can  be  calculated  as 

y*  =  AQ  +  AjX  (11) 

2  2 

with  an  error,  e  ,  calculated  as  e  =  (y  -  y*) 

61.  Using  this  procedure,  a  regression  equation  is  defined  which,  when 
applied,  has  a  zero  mean  error  and  a  minimum  error  variance.  The  same  proce¬ 
dure  outlined  above  is  used  for  linear  equations  of  more  than  one  independent 
variable  or  for  nonlinear  functional  equations.  A  modification  is  required, 
however,  to  accommodate  additional  variables  or  to  introduce  the  nonlinear 
functional  form. 


Regression  Analysis  for  Estimation  o f_ the  Joint 
Roughness  Coefficient 

62.  Plots  of  JRC  versus  fractal  dimension  and  JRC  versus  variogram 
nugget  are  shown  in  Figure  43.  A  linear  correlation  is  observed  between  the 
JRC  and  each  of  the  independent  variables.  Based  on  this  preliminary  visual 
analysis,  three  forms  of  a  regression  equation  are  proposed: 


a. 

JRC  =  A 

+ 

A,D. 

0 

1 

b. 

JRC  =  An 

+ 

A.  CO. 

0 

1 

c . 

JRC  =  AQ 

+ 

AjD  +  A2CO. 

where 

D  =  fractal  dimension 
CO  =  variogram  nugget 

Each  of  these  functional  forms  was  tested  to  determine  the  optimal  regression 
equation. 

63.  The  data  used  to  define  the  coefficients  of  these  equations  are 
presented  in  Table  9.  The  fictitious,  ideal  cross  section  was  used  to  include 
the  ideal  case  where  JRC  is  equal  to  zero.  Each  equation  was  evaluated  for 
accuracy  relative  to  the  mean  square  error  resulting  from  its  application. 
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Figure  43.  Two  regression  plots:  JRC  versus  fractal  dimension  and 

JRC  versus  variogram  nugget 
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64.  Using  the  method  of  least  squares,  the  following  regression  equa¬ 
tions  were  derived: 

a.  For  the  fractal  dimension  only: 

JRC  =  -1022.55  +  (1023. 92)D  (19) 

b.  For  the  variogram  nugget  only: 

JRC  =  2.042  +  (818 . 53)CO  (20) 

£.  For  the  fractal  dimension  and  variogram  nugget: 

JRC  =  -758.41  +  (759 . 93)D  +  (215.05)CO  (21) 

For  each  of  these  equations,  D  and  CO  must  be  represented  by  at  least  four 
digit  precision  after  the  decimal  point.  Table  10  lists  the  results  of  the 
application  of  each  equation  for  the  estimation  of  JRC.  The  most  accurate 
equation  is  the  one  based  on  the  two  independent  variables. 

65.  It  is  acknowledged  that  the  accuracy  of  the  three  equations  is 
similar.  Moreover,  these  equations  are  empirical  and  are  forwarded  only  as  a 
procedure  for  the  estimation  of  JRC.  It  is  difficult  for  users  unfamiliar 
with  variograms  to  interpret  the  nugget  value.  It  is  therefore  recommended  to 
use  the  regression  equation  based  solely  on  the  fractal  dimension  to  estimate 
JRC.  This  equation,  in  review,  is 

JRC  =  -1022.55  +  (1023. 92)D  (19) 

This  equation  has  an  error  standard  deviation  equal  to  1.24  JRC  values.  When 
using  this  equation,  it  is  essential  that  the  fractal  dimension,  D  ,  be 
represented  by  at  least  four  digit  precision  after  the  decimal  point,  even 
though  the  coefficients  have  a  precision  less  than  this.  For  rock  surfaces, 
the  fractal  dimension  has  a  range  from  1  to  about  1.03;  hence,  these  values 
are  small  and  the  precision  is  therefore  necessary. 

66.  A  rough  approximation  of  equation  (19)  is  developed  as 


JRC  =  1000(D  -  1) 


(22) 


Table  10 

Application  of  Regression  Equations  for  the  Estimation  of  JRC 


Rock  Surface 

JRC 

JRC1* 

JRC2 

3 

JRC 

Fictitious  Ideal 

0 

1.37 

2.04 

1.52 

Backus  Notch  No.  1  (bedding) 

3 

2.50 

2.78 

2.55 

Backus  Creek  Notch  No.  1  (bedding) 

2 

2.29 

2.70 

2.38 

DS  +  122  (C  joint) 

4 

4.54 

5.23 

4.71 

Island  Notch  (joint) 

1 

1.47 

2.12 

1.62 

DS  +  122  (minor  joint) 

4 

4.85 

4.17 

4.66 

DS  +  122  (joint  No.  2) 

5 

4.95 

4.25 

4.76 

914  Rib  (joint) 

3 

4.75 

4.33 

4.63 

Dunn  Creek  Notch  (joint) 

7 

6.18 

4.91 

5.84 

930  (joint) 

8 

5.98 

6.13 

6.01 

Island  Notch  (joint/bed) 

5 

4.54 

4.09 

4.41 

Island  Notch  (bedding) 

2 

2.80 

3.35 

2.93 

Wolf  Creek  Jet.  (joint) 

8 

7.10 

5.32 

6.64 

Dunn  Creek  Notch  (bedding) 

13 

10.18 

13.01 

10.94 

Old  Notch  No.  1  (joint) 

22 

23.90 

22.83 

23.70 

Backus  Notch  No.  1  (joint) 

6 

5.67 

5.73 

5.68 

Mean  Square  Error 

1.54 

1.77 

1.39 

Mean  Standard  Deviation 

1.24 

1.33 

1.18 

*  Superscripts  on  the  JRC  headings  represent 

the  equation 

used  for 

estima- 

tion.  Superscript  1  represents  the 

equation 

based  on  the 

fractal 

dimension 

alone;  superscript  2  represents  the 

equation  based  on  the 

variogram  nugget 

alone;  and  superscript  3  represents 

the  equation  based  on 

the  fractal  dimen- 

sion  and  variogram  nugget. 

and  this  equation  is  more  conservative 

s  than  equation  (19); 

that  is, 

a  lower 

estimate  of  JRC  results  from  equation 

(22).  Moreover,  this 

equation 

,  is  useful 

for  the  calculation  of  other  rock  surface  parameters,  such  as  shear  stress. 
For  example,  equation  (9)  is  recalled: 


JRC  = 


arctant  (x/cr  )  -  $ 

N  b 

Log10(JCS)/aN) 


(9) 


If  the  shear  stress,  t  ,  is  to  be  calculated  based  on  a  known  JRC  value,  then 


tan-1(T/aN)  =  JRC(log10( JCS/oN) )  + 


(23) 
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or 

t  =  oNtan(JRC(log10(JCS/aN>)  +  <J>b)  (24) 

Incorporating  the  fractal  dimension,  D  ,  into  the  equation  is  most  easily 
done  using  equation  (22): 

t  =  aNtan(I00(D  -  1)  (log1()(  JCS/cy  )  +  4>b)  (25) 

Equation  (25)  is  an  expression  for  the  calculation  of  shear  stress  as  a  func¬ 
tion  of  the  fractal  dimension. 

67.  An  additional  study  was  made  to  investigate  the  accuracy  of  equa¬ 
tions  (19)  and  (22).  An  enlargement  of  Figure  41  was  made,  and  a  fractal 
dimension  was  computed  for  each  of  the  ten  profiles.  The  results  of  these 
calculations  are  presented  in  Table  11  with  a  plot  of  JRC  versus  fractal  dimen¬ 
sion  shown  in  Figure  44.  It  is  interesting  that  through  the  JRC  range  of  14, 
the  plot  of  JRC  versus  fractal  dimension  is  steeper  in  Figure  44  than  what  is 
shown  by  Figure  43;  i.e.,  based  on  Figure  41,  an  increase  in  the  fractal 
dimension  should  indicate  a  larger  increase  in  the  JRC  value  than  what  is 
predicted  by  equations  (19)  and  (22).  For  the  JRC  range  greater  than  14,  how¬ 
ever,  Figure  44  shows  a  distinct  change  in  the  slope  of  the  plot  of  JRC  versus 
fractal  dimension.  Several  explanations  are  possible.  The  resolution  of  the 
enlargement  (a  photographic  enlargement)  was  poor,  and  this  might  account  for 
the  change.  Barton  and  Choubey  (1977)  seemed  to  assign  high  JRC  values  to 
profiles  which  had  a  sinusoidal  shape  but  were  otherwise  smooth  (i.e.,  similar 
to  a  sine  curve).  A  low  fractal  dimension  would  result  for  such  a  curve. 

Some  of  the  profiles  for  higher  JRC  values  appear  to  have  such  a  shape  (Fig¬ 
ure  41).  It  is  observed  from  Table  12  that  the  application  of  equations  (19) 
and  (22)  to  the  profiles  of  Figure  41  result  in  lower  estimates  of  JRC  than 
what  is  shown  in  Table  11.  Thus,  equations  (19)  and  (22)  are  useful  for  a 
conservative  analysis  of  joint  surface  roughness. 
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Table  11 

Fractal  Dimension  Calculations  for  Profiles  of  Figure  41 


JRC 

Number 

of  Segments, 

N 

Fractal 

Range 

y  =  1  cm 

y  =  2  cm 

y  =  3  cm 

y  =  4  cm 

y  =  5  cm 

Dimension 

0-2 

173.267 

86.589 

57.708 

43.277 

34.621 

1.000603 

2-4 

178.600 

89.242 

59.465 

44.570 

35.656 

1.001199 

4-6 

178.408 

89.047 

59.332 

44.450 

35.552 

1.002287 

6-8 

176.569 

88.117 

58.691 

43.948 

35.138 

1.003124 

8-10 

176.468 

88.036 

58.616 

43.893 

35.073 

1.003822 

10-12 

177.890 

88.701 

59.028 

44.207 

35.301 

1.004689 

12-14 

172.990 

86.251 

57.401 

42.963 

34.317 

1.005596 

14-16 

184.614 

92.123 

61.238 

45.788 

36.560 

1.006203 

16-18 

182.133 

90.697 

60.244 

44.942 

35.874 

1.009647 

18-20 

184.560 

91.550 

60.686 

45.321 

36.205 

1.012453 

Table  12 

Predicted  JRC  Values  for  Profiles  of  Figure  41 


JRC  Range 

Fractal 

Predicted 

JRC  Values 

(actual) 

Dimension 

Eq.  (19) 

Eq.  (22) 

0-2 

1.000603 

1.99 

0.60 

2-4 

1.001199 

2.60 

1.20 

4-6 

1.002287 

3.71 

2.29 

6-8 

1.003124 

4.57 

3.12 

8-10 

1.003822 

5.28 

3.82 

10-12 

1.004689 

6.17 

4.69 

12-14 

1.005596 

7.10 

5.60 

14-16 

1.006203 

7.72 

6.20 

16-18 

1.009647 

11.25 

9.65 

18-20 

1.012453 

14.12 

12.45 
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PART  VII:  CONCLUSIONS 


Implication  of  Results 

68.  Applications  of  both  the  fractal  dimension  and  the  variogram  were 
novel  for  the  analysis  of  the  texture  of  a  joint  surface.  The  objective  of 
this  study  was  to  evaluate  each  of  these  techniques  for  the  description  of 
rock  surfaces.  As  was  reviewed  in  this  report,  both  of  these  techniques  are 
valuable  for  rock  mass  characterization. 

69.  The  fractal  dimension  accurately  measured  the  roughness  of  discrete 
cross  sections  across  the  DS  +  122  joint  surface.  Regions  which  appeared  to 
be  smooth  on  the  contour  map  were  found  to  have  a  lower  fractal  dimension  com¬ 
pared  with  regions  of  greater  complexity.  Moreover,  the  value  of  the  fractal 
dimension  is  directly  proportional  to  surface  roughness.  The  greater  the 
difference,  D  -  1  ,  the  rougher  is  the  surface.  The  fractal  dimension  offers 
a  unique  approach  to  the  quantitative  description  of  surface  roughness,  which 
has  otherwise  been  qualitatively  described.  These  results  were  also  verified 
for  the  string  line  cross  sections. 

70.  As  a  complement  to  the  fractal  dimension,  variograms  for  selected 
cross  sections  demonstrated  that  elevation  values  were  spatially  correlated. 
The  variance  in  elevation  values  along  each  cross  section  presented  in  this 
report  is  a  function  of  the  asperities  along  the  joint  surface.  Variograms 
were  able  to  show  that  asperities  are  spatially  correlated  and  have  a  finite 
spatial  continuity.  Moreover,  at  least  for  the  example  of  Figure  25,  the 
variogram  is  useful  for  the  investigation  of  asperity  magnitude  differences 
across  a  joint  surface.  Occasionally,  the  variogram  yields  ambiguous  results 
concerning  asperity  magnitude  and  spatial  structure.  This  ambiguity  was 
evident  for  some  of  the  string  line  cross  sections. 

Problems  Encountered  in  this  Research 


71.  The  only  problems  which  occurred  during  this  research  resulted  from 
the  application  of  existing  techniques  to  a  novel  problem.  The  concept  of  the 
fractal  dimension  is  a  recent  development.  The  literature,  except  for  publi¬ 
cations  by  Mandelbrot,  only  sparsely  contains  publications  whic!  discuss  the 
fractal  dimension.  It  is  therefore  difficult,  at  present,  to  check  the 
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correctness  of  a  research  procedure  which  involves  the  calculation  of  a 
fractal  dimension.  Such  verification  was  provided  in  this  report  using  a  con¬ 
trol  calculation. 

72.  With  respect  to  the  variogram,  geostatistical  techniques  are  based 
on  the  principle  oft  he  Intrinsic  Hypothesis.  This  hypothesis  simply  requires 
that  the  local  mean  of  the  data  in  one  location  be  identical  to  the  local  mean 
in  a  different  location.  This  is  an  hypothesis  required  by  the  geostatistical 
techniques.  For  this  study,  the  elevation  data  for  some  of  the  cross  sections 
conformed  to  the  hypothesis.  Elevation  values  for  other  cross  sections,  how¬ 
ever,  did  not  conform  well  to  this  hypothesis.  In  these  latter  cases,  it  was 
necessary  to  filter  the  data  to  allow  a  conformity  with  the  Intrinsic  Hypothe¬ 
sis.  This  was  not  necessary  for  the  string  line  data.  These  data  were  sam¬ 
pled  relative  to  a  string  stretched  above  and  parallel  to  the  rock  surface. 
This  method  of  data  collection  yields  stationary  data. 

Solutions  to  Research  Problems 


73.  For  the  calculation  of  the  fractal  dimension,  some  experimentation 
is  required.  The  better  approach  to  this  calculation  involves  four  steps: 

(a)  select  a  finite,  small  segment  length,  y  ;  (b)  count  the  number  of 
times,  N  ,  the  segment  length,  y  ,  can  be  laid,  end  to  end,  along  an  irregu¬ 
lar  line;  at  the  end  of  this  irregular  line,  add  the  remainder,  k  ,  norma¬ 
lized  by  y  to  the  counts,  N  (i.e.,  N  =  N  +  k/y);  (c)  repeat  steps  a  and  b 
for  different  segment  lengths,  y  ;  (d)  plot  log^(N)  versus  log^(y);  the 
negative  of  the  slope  of  this  line  is  the  fractal  dimension.  As  an  alterna¬ 
tive,  the  fractal  dimension  can  be  calculated  using  the  computer  program 
listed  in  Appendix  E. 


Research  Contribution 


74.  Once  the  fractal  dimension  has  been  calculated,  a  JRC  value  can  be 
estimated  using  the  fractal  dimension.  This  study  found  that  JRC  and  the 
fractal  dimension  are  strongly  correlated.  This  is  a  contribution  to  the 
state-of-the-art  in  rock  mechanics  because  a  procedure  is  now  available  to 
estimate,  or  calculate,  the  JRC  rather  than  guessing  its  value  from  a  figure, 
such  as  Figure  41.  The  relationship  between  JRC  and  the  fractal  dimension  is 
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JRC  -  -1022.55  +  (1023. 92)D 


(19) 


with  D  ,  the  fractal  dimension,  having  a  precision  of  at  least  four  digits 
after  the  decimal  point.  Such  precision  is  easily  obtained  from  the  computer 
program  listed  in  Appendix  E. 
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APPENDIX  A:  THE  CONCEPT  OF  THE  FRACTAL  DIMENSION 


Introduction 


1.  Geometry  is  the  branch  of  mathematics  that  addresses  the  measurement, 
properties,  and  relationships  of  points,  lines,  angles,  surfaces  and  solids. 

In  a  broad  context,  this  branch  of  mathematics  studies  the  properties  of  given 
elements  which  remain  invariant  under  specified  transformations.  Geometry 
also  represents  an  arrangement  of  objects  or  parts  which  suggest  a  particular 
shape  or  figure. 

2.  In  nature,  shapes  are  geometrically  unique.  The  shapes  of  the  Amazon 
drainage  system,  a  butterfly's  wings,  a  blood  cell,  or  a  cubic  crystal  of 
galena  are  unique.  Moreover,  these  shapes  are  often  difficult  to  describe. 

The  crystal  of  galena  is  geometrically  familiar:  a  cube.  A  blood  cell, 
although  more  complex  in  shape  than  a  cube,  is  described  well  as  a  thin,  hol¬ 
lowed  cylinder  or  disk.  In  contrast,  there  is  no  easy  geometrical  way  to 
describe  the  shape  of  the  Amazon  drainage  system. 

3.  From  a  topological  viewpoint,  the  shape  of  the  Amazon  River  is  a 
curve.  The  same  is  true  of  any  coastline.  A  curve  is  a  line,  a  one  dimen¬ 
sional  entity.  Hence,  the  topological  dimension  of  a  river  system  or  a  coast¬ 
line  is  one.  This  characterization  is,  however,  disturbing.  The  Amazon 
drainage  system  does  not  look  like  a  familiar  curve  or  line.  It  looks  instead 
like  a  complex  system  of  smaller  curves  and  lines,  randomly  oriented  into  a 
dendritic  collage.  As  a  result  of  this  example,  a  topological  description  of 
a  river  system  is  incomplete. 


Fractal  Dimension 


4.  Rather  than  attempting  to  describe  the  shape  of  the  Amazon  River 
system,  or  any  other  natural  and  irregular  form,  an  effort  can  be  made  to 
classify  this  irregularity.  This  effort  can  be  achieved  through  the  attempt 
to  measure  the  length  of  the  irregular  form,  in  this  case,  the  length  of  the 
Amazon  River. 

5.  The  Amazon  River  is  shown  to  scale  in  Figure  Al.  The  length  of  this 
stream  can  be  approximated  by  selecting  a  particular  ruler  length,  y  ,  then 


Al 


counting  the  number  of  times,  N  ,  the  length,  y  ,  can  be  laid  along  the 

stream.  The  length  of  the  stream,  L  ,  is  approximated  as  Ny  . 

6.  If  the  ruler  length,  y  ,  is  progressively  decreased  and  allowed  to 
converge  to  zero,  the  length,  L  ,  of  the  stream  is  found  to  increase  without 
bound  (Mandelbrot  1967)*.  This  statement  seems  to  be  incorrect  because  it  is 
thought  that  L  converges  to  some  finite  value  as  y  converges  to  zero. 

This  is  not  the  case,  however,  as  is  graphically  demonstrated  in  Figure  A2. 
Imagine  if  one  were  to  measure  the  length  of  the  Amazon  in  situ,  not  from  a 
map  trace.  If  y  is  1  km,  one  result  for  L  would  be  obtained.  As  y  is 

decreased,  say  to  1  m,  a  different  result  for  L  would  be  obtained  because 

features,  such  as  embankment  undulations,  become  important  at  the  scale  of  1  m 
that  were  overlooked  at  the  scale  of  1  km.  As  y  is  further  decreased,  say 
to  1  mm,  individual  rock  fragments  in  the  embankment  of  the  river  become 
important  and  add  to  the  length,  L  .  This  describes  the  unbounded  increase 
in  the  length,  L  . 

7.  Because  of  the  unbounded  increase  of  the  global  length,  L  ,  this 
parameter  is  not  satisfactory  for  the  description  of  natural  shapes.  The 
relationship  between  the  ruler  length,  y  ,  and  the  global  length,  L  ,  is 
valuable,  however,  for  the  description  of  natural  shapes.  This  relationship 
was  given  before  as 

L  =  Ny  (Al) 

This  relationship  can  be  expressed  in  a  slightly  different  fashion  as 

L  =  Ny°  (A2) 

where  D  is  called  the  Hausdorf f-Besicovitch  dimension  (Mandelbrot  1982).  In 
equation  (Al),  this  dimension  is  equal  to  1. 

8.  For  the  Amazon  River,  each  of  the  ruler  lengths,  y  ,  (1  km,  1  m, 

1  mm),  is  associated  with  a  particular  number  of  steps,  N  ,  used  to  approxi¬ 
mate  the  global  length,  L  .  Figures  A3  through  A5  show  the  results  for  the 
Amazon  River  less  its  tributaries  with  ruler  lengths  of  160  km,  80  km  and 

*  References  cited  in  this  appendix  are  included  in  the  references  at  the  end 
of  the  main  text. 
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Figure  A3.  Amazon  River  trace,  showing  the  measurement  of  its  length 

using  160  km  step  sizes  (y) 


A6 


Figure  A4 .  Amazon  River  trace,  showing  the  measurement  of  its  length 
using  80  km  step  sizes  (y) .  Here,  20  +  K  steps  are  shown,  whe^e  k 


A7 


Figure  A5.  Amazon  River  trace,  showing  the  measurement  of  its  length 
using  25  km  step  sizes  (y)  .  Here,  66  +  k  steps  are  shown,  where  K 


25  km.  As  y  decreases,  N  increases  without  bound;  hence,  L  increases 
without  bound. 

9.  By  examining  equation  (A2) ,  however,  and  rearranging  terms,  the 
following  relationship  is  found: 


Ly  D  =  N  (A3) 

and  there  is  a  constant  relationship  between  y  and  N  .  If  the  base-ten  log 
of  N  is  plotted  against  the  base-ten  log  of  y  ,  the  slope  of  this  plot  is 
-D  .  This  relationship  is  expressed  by  equation  (A3).  Moreover,  the  slope  of 
this  plot  is  a  straight  line;  hence,  it  is  constant.  Therefore,  despite  the 
divergent  behavior  of  the  global  length,  L  ,  the  Hausdorf f-Besicovitch  dimen¬ 
sion  is  constant. 

10.  The  results  from  Figures  A3  through  A5  are  plotted  in  Figure  A6. 
From  this  plot,  it  is  evident  that  the  relationship  between  log^(N)  and 
log^(y)  is  constant.  The  slope  of  this  line  is  -1.02.  Therefore,  the 
Hausdorf f-Besicovitch  dimension  is  1.02. 

11.  If  the  Amazon  River  were  perfectly  straight,  as  shown  in  Figure  A7 , 
a  different  result  would  be  obtained  for  the  Hausdorf f-Besicovitch  dimension. 
The  log^(N)  versus  log^Cy)  plot  for  this  idealized  river  is  shown  in  Fig¬ 
ure  A8.  The  slope  of  this  plot  is  exactly  -1  ,  thus,  the  Hausdorff- 
Besicovitch  dimension  is  1 . 

12.  The  hypothetical  stream  of  Figure  A7  is  represented  perfectly  by 
the  topological  line.  For  this  stream,  D  =  1;  therefore,  equation  (Al)  accu¬ 
rately  measures  the  length.  As  was  stated  before,  the  topological  dimension 
of  any  river  is  1.  For  the  stream  of  Figure  A7,  the  Hausdorf f-Besicovitch 
dimension  is  exactly  equivalent  to  the  topological  dimension. 

13.  This  is  not  the  case  for  the  Amazon  River.  Its  Hausdorff- 
Besicovitch  dimension  was  found  to  be  1.02.  This  value  is  greater  than  1.  In 
nature,  any  curve,  surface,  or  volume  whose  Hausdorf f-Besicovitch  dimension 
exceeds  its  topological  dimension  is  a  fractal.  Furthermore,  the  Hausdorf f- 
Besicovitch  dimension  is  known  as  the  fractal  dimension  (Mandelbrot  1982). 

14.  The  fractal  dimension  provides  the  capability  for  the  description 
of  irregular,  natural  shapes.  Whereas  the  fractal  dimension  is  simply  a 
numerical  value  and  does  not  describe  a  shape  such  as  a  line,  it  does  describe 
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Figure  A6.  Log^(N)  versus  Log^(y)  plot  for  the  Amazon  River. 
The  fractal  dimension,  D  ,  is  1.02 
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the  complexity  or  irregularity  of  a  curve,  a  surface,  or  a  volume.  This 
offers  considerable  capability  for  the  description  of  natural  surfaces. 

Example  Calculations 

15.  An  example  calculation  of  the  fractal  dimension  has  been  provided 
in  Figures  A3  through  A6.  The  concept  of  the  fractal  dimension  has  also  been 
applied  previously  for  the  description  of  coastlines  (Mandelbrot  1967). 

16.  Two  coastlines  are  chosen  as  examples.  Figure  A9  shows  the  west 
coast  of  Great  Britain;  Figure  A10  shows  the  coast  of  southern  California. 

The  southern  California  coast  qualitatively  appears  to  be  smoother  than  the 
west  coast  of  Great  Britain.  Therefore,  it  is  expected  that  the  fractal 
dimension  of  the  southern  California  coast  will  be  less  than  the  fractal 
dimension  for  the  west  coast  of  Great  Britain. 

17.  Using  the  scales  plotted  on  each  map,  the  selection  of  the  ruler 
length,  y  ,  is  as  shown  in  Tables  A1  and  A2.  The  technique  for  determining 
N  for  each  y  is  the  same  as  was  reviewed  in  Figures  A3  through  A5.  The 
only  difference  in  the  procedure  involved  the  use  of  different  ruler  lengths, 

y  • 

18.  From  Tables  A1  and  A2 ,  the  log^CN)  versus  log^(y)  plots  for  each 
coastline  are  shown  in  Figures  All  and  A12.  From  these  plots,  the  fractal 
dimension  for  the  west  coast  of  Great  Britain  is  1.25.  The  fractal  dimension 
of  the  southern  California  coast  is  1.02.  Hence,  the  fractal  dimension  veri¬ 
fied  the  intuition  that  the  southern  California  coast  is  smoother  than  the 
west  coast  of  Great  Britain. 


Mandelbrot  Set 


19.  The  Mandelbrot  Set  is  a  recent  discovery  in  mathematics.  This  set 
consists  of  a  region  at  the  contact  between  the  plane  of  real  numbers  and  the 
plane  of  complex  numbers.  Within  this  region,  numerical  vectors  defined  by 
the  real  and  complex  planes  remain  of  finite  numerical  value  following  an 
infinite  number  of  iterations  in  which  the  value  of  each  vector  is  consis¬ 
tently  squared. 
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Figure  A9.  Outline  of  the  coast  of  Great  Britain 


Figure  A10.  Outline  of  the  coast  of  southern  California 


Table  A1 

Fractal  Dimension*  Calculation,  West  Coast  of  Great  Britain 


Ruler  Length 
y  (km) 

logio(y) 

Segment  Counts,  N 

logioW 

125 

2.097 

7.57 

0.879 

50 

1.699 

24.25 

1.385 

25 

1.398 

56.10 

1.749 

*  Fractal  Dimension  is  1.25. 


Table  A2 

Fractal  Dimension*  Calculation,  Southern  California  Coast 


Ruler  Length 
y  (km) 

M 

O 

OQ 

O 

/*\ 

Segment  Counts,  N 

log10 W 

100 

2.000 

7.600 

0.880 

40 

1.600 

19.330 

1.290 

10 

1.000 

80.000 

1.903 

*  Fractal  Dimension  is  1.02. 
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Figure  All.  Logjg(N)  versus  Log^(y)  plot  for  the  western  coast  of 
Great  Britain.  The  fractal  dimension  is  1.25 
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Figure  A12.  Log^(N)  versus  Log^Cy)  plot  for  the  southern  California 
coast.  The  fractal  dimension  is  1.02 
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20.  An  example  of  the  fringe  surrounding  the  Mandelbrot  Set  is  shown  in 
Figure  A13.  The  fringe  (shaded  region)  is  not  part  of  the  Mandelbrot  Set; 
this  set  resides  inside  this  fringe.  The  IBM-PC  compatible,  BASIC  language 
program  used  to  create  this  figure  (Dewdney  1985)  is  listed  in  Figure  A14. 

21.  The  Mandelbrot  Set  is  a  fractal.  This  set  is  not  the  result  of  a 
contrived  mathematical  function  used  to  create  an  interesting  graphical  image. 
This  set  is,  instead,  a  fascinating  property  of  real  and  complex  numbers. 

This  set  provides  one  more  example  of  the  fractal  geometry  of  nature. 


-2.0  '.75  +.50 


REAL  PART 

Figure  A13.  Mandelbrot  Set.  The  border  shows  the  range  of 
real  and  imaginary  numbers  used  to  create  this  object. 
(ACORNER  =  -2.0,  BCORNER  =  1.25,  SIDE  =2.5) 
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5  KEY  OFF:  SCREEN  1:  COLOR  0.0 
10  INPUT  “ENTER  THE  VALUE:  ACORNER" ; 
20  INPUT  “ENTER  THE  VALUE:  BCORNER" ; 

50  INPUT  "ENTER  THE  VALUE:  SIDE 

55  GAP  =  SIDE  /  150.0 

56  CLS 

60  FOR  I  =  1  TO  150 

70  FOR  J  =  1  TO  150 

90  C  =  <J  *  GAP  +  ACORNER) 

95  E  =  (BCORNER  -  I  *  GAP) 

96  R  -  C 

97  W  -  E 

100  I COUNT  =  0 


101 

U  =  <R"2  - 

W  2 

+  C) 

1  02 

V  =  (2*R*W 

+  E) 

no 

Z  =  SQR  ( U  "2 

+  V 

-2) 

1 20 

I  COUNT  =  I COUNT 

+  1 

1 30 

IF  Z  >  2  THEN  150 

135 

IF  I COUNT  > 

90 

THEN 

150 

136 

R  =  U 

137 

W  =  V 

138 

GO  TO  101 

150 

L  =  0 

151 

IF  I COUNT  < 

10 

THEN 

164 

152 

IF  I COUNT  > 

16 

THEN 

155 

153 

L  =  1 

154 

GO  TO  164 

155 

IF  I COUNT  > 

•-}  ~T 

THEN 

GO 

ifJ 

156 

L  =  2 

157 

GO  TO  164 

158 

IF  I COUNT  > 

30 

THEN 

164 

159 

L  =  3 

164 

F'SET  (J.I). 

L 

1 70 

NEXT  J 

175 

NEXT  I 

130  STOP 


Figure  A14.  IBM  compatible  BASIC  language  program 
the  Mandelbrot  Set  of  Figure  A13 


ACORNER 

BCORNER 

SIDE 


used  to  create 
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APPENDIX  B:  THE  THEORY  OF  REGIONALIZED  VARIABLES, 
A  GEO STATISTICAL  TECHNIQUE 


1.  The  theory  of  regionalized  variables  was  developed  by  Georges 
Matheron  in  the  late  1950's.  Matheron  demonstrated  that  spatially  dependent 
variables  can  be  estimated  on  the  basis  of  this  spatial  structure  and  known 
samples  (Matheron  1963)*.  This  estimation  is  one  aspect  of  geostatistics,  a 
concept  concerned  with  describing  the  distribution,  in  space,  of  geologic 
phenomena . 

2.  A  random  variable  distributed  in  space  is  said  to  be  regionalized. 
These  variables,  because  of  their  spatial  aspect,  possess  both  random  and 
structured  components.  On  a  local  scale,  regionalized  variables  can  behave 
randomly  or  erratically.  Two  regionalized  variables  separated  by  a  distance, 
h  ,  however,  are  not  independent,  but  are  related  by  a  structured  aspect 
dependent  upon  h  .  Usually,  as  the  length  of  h  increases,  the  similarity 
between  two  regionalized  variables  decreases. 

3.  At  first  glance,  a  regionalized  variable  appears  to  be  a  contradic¬ 
tion.  In  one  sense,  it  is  a  random  variable  which  locally  may  have  no  rela¬ 
tion  to  surrounding  variables.  On  the  other  hand,  there  is  a  structured 
aspect  to  a  regionalized  variable  which  depends  on  the  distance  separating  the 
variables.  Both  of  these  characteristics  can,  nevertheless,  be  described  by  a 
random  function  of  which  each  regionalized  variable  is  a  single  realization. 

By  incorporating  both  the  random  and  structured  aspects  of  a  regionalized 
variable  in  a  single  function,  spatial  variability  can  be  accommodated  on  the 
basis  of  the  spatial  structure  shown  by  these  variables. 

Variogram 

4.  One  way  to  examine  the  spatial  structure  of  a  regionalized  variable 
is  to  analytically  relate  the  change  in  samples,  or  measurements,  of  the  vari¬ 
able  as  a  function  of  distance  separating  the  samples.  In  general,  if  the 
average  difference  between  samples  increases  as  their  distance  of  separation 
increases,  a  spatial  structure  exists,  and  the  variable  is  regionalized. 


*  References  cited  in  this  appendix  are  included  in  the  references  at  the  end 
of  the  main  text. 
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5.  The  function  which  defines  the  spatial  correlation,  or  structure,  of 
a  regionalized  variable  is  the  variogram.  This  function  is  denoted  as  y  and 
is  defined  as 


N 

Y(h)  =  2N  ^ 
i=l 

where  N  is  the  total  number  of  data  pairs  separated  by  a  distance,  h  .  The 
variogram  is  simply  one  half  the  average  square  of  the  difference  between  sam¬ 
ples,  Z(X^)  ,  separated  by  a  distance,  h  .  If  a  spatial  relationship 
exists,  the  value  of  y(h)  increases  as  the  separation  distance,  h  , 
increases.  This  implies  that  samples  located  close  in  space  are  more  similar 
in  value  than  those  separated  by  a  considerable  distance. 

6.  For  most  geostatistical  applications,  the  variogram  has  a  spherical 
shape  (Journel  and  Huijbregts,  1978).  This  function  can  be  modeled  by  the 
following  equation: 


z(xi) 


Cri 


Nugget  +  C  ^ 


P.5h3 

R3 


,  h  <  R 


Y  (h)  = 

Nugget  +  C,  h  £  R 


(B2) 


Each  component  of  this  equation  is  illustrated  in  Figure  Bl.  The  "nugget" 
value  indicates  the  amount  of  white,  or  random,  noise  present  in  a  set  of 
data.  This  value  is  the  intercept  of  the  variogram  with  the  abscissa  axis. 
The  value,  R  ,  is  the  distance  at  which  the  variogram  stops  increasing  and 
becomes  flat;  this  value  is  known  as  the  range.  The  quantity,  C  ,  is 
obtained  by  subtracting  the  nugget  from  the  sill.  The  sill  is  the  variogram 
value  at  the  distance,  h  =  R  ,  and  is  usually  equal  to  the  sample  variance. 


Linear  Estimation  of  the  Regionalized  Variables:  Kriging 

7.  Once  the  spatial  structure  of  a  regionalized  variable  has  been 
determined  through  computation  of  the  variogram,  the  spatial  structure  can  be 
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Figure  Bl.  Example  of  a  spherical  variogram  showing 
the  nugget,  range,  sill  and  shape 
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used  to  estimate  the  value  of  a  variable  at  an  unsampled  location.  This  esti¬ 
mation,  or  interpolation  process,  is  known  as  kriging  (Matheron,  1963). 

8.  The  estimate,  Z*(X^)  ,  of  a  regionalized  variable  at  a  location, 

XQ  ,  is  given  by 


N 

z*(v  xiz(xi>  (b3) 

i=l 


where  N  is  the  number  of  closest  surrounding  locations  used  for  estimation. 

9.  In  equation  (B3)  ,  A  is  the  vector  of  weights  applied  to  each  of 
the  closest  surrounding  locations.  Generally,  the  largest  weights  are 
assigned  to  the  closest  locations.  Kriging  is  an  unbiased  estimator  and 
results  in  a  minimum  variance  of  the  error.  These  objectives  are  controlled 
by  the  weights,  A  .  For  the  unbiased  condition,  the  sum  of  the  weights  is 
unity: 


N 

\  =  1  (B4) 

i=l 

where  N  is  the  same  as  for  equation  (B3) .  To  compute  the  vector,  A  ,  the 
following  procedure  is  used: 


E 


-  V 


(B5) 


Where  p  is  an  unknown,  solved  along 
plier  used  to  achieve  a  minimum  error 
are  respectively  the  intersample 
obtained  by  solving 


with  A  ,  and  is  the  Lagrangian  multi¬ 
variance.  In  equation  (B5) ,  o  and 
and  point-sample  covariances.  These  are 


o(h)  =  Sill  -  y (h) 


(B6) 


For  the  intersample 
rating  each  of  the 


covariance  matrix, 
N  known  locations 


o^  ,  the  distance, 
used  for  estimation. 


h  ,  is  that  sepa- 
For  the 
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point-sample  covariance  vector,  o  ,  •  the  distance,  h  ,  is  that  between  the 

oi 

location,  ,  and  each  of  the  N  known  locations  used  for  estimation.  In 
equation  (B6) ,  y  is  obtained  from  equation  (B2). 

10.  Once  the  weighting  vector,  X  ,  is  known,  the  kriging  variance,  which 
is  the  variance  of  the  estimation  error,  is  computed  as 

N 

Krig  Var  =  Sill  -  ^  Vio"  p  (B7) 

i=l 


where  the  terms  sill,  X  ,  cr  ^  ,  and  y  are  defined  as  above.  In  a  practi¬ 
cal  sense,  the  kriging  variance  is  analogous  to  the  mean  square  error  of  the 
estimate. 


Intrinsic  Hypothesis 

11.  This  hypothesis  is  also  known  as  the  concept  of  data  stationarity 

(Journel  and  Huijbregts  1978) .  The  Intrinsic  Hypothesis  holds  that  the  mean 
of  the  data  in  one  location,  ,  is  equivalent  to  the  mean  at  another  loca¬ 

tion,  X^  .  This  requires  that  data  values  be  more  or  less  equal  throughout  a 
space,  except  for  their  natural  variation  (variance). 

12.  Kriging  is  a  technique  which  expects  each  estimate  to  be  the  mean  of 
the  data.  This  explains  why  kriging  is  most  accurate  where  the  Intrinsic 
Hypothesis  holds.  The  natural  variability  in  the  data  is  captured  by  the 
variogram.  Thus,  kriging  is  really  attempting  to  model  data  fluctuations. 

13.  The  Intrinsic  Hypothesis  cannot  hold  everywhere  in  space;  otherwise, 
a  natural  phenomenon  would  be  ubiquitous  throughout  the  universe.  Eventually, 
a  natural  phenomenon  diminishes  and  becomes  non-existent.  An  ore  body  gradu¬ 
ally  grades  to  barren  country  rock;  earthquake  ground  motion  gradually  dimin¬ 
ishes  away  from  an  epicenter;  and  the  concentration  of  oxygen  in  the  earth’s 
atmosphere  gradually  lessens  as  one  ascends  through  the  atmosphere.  These  are 
all  examples  of  nonstationary  data  behavior.  That  is,  the  mean  is  not  con¬ 
stant,  but  gradually  decreases.  For  these  examples,  the  Intrinsic  Hypothesis 
does  not  hold. 

14.  Despite  the  examples  of  non-stationarv  behavior,  the  data  will  still 
show  natural  fluctuations.  By  invoking  the  concept  of  quasi-stationary 
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behavior  (Journel  and  Huijbregts  1978),  we  can  restrict  the  size  of  the  space 
over  which  we  require  the  Intrinsic  Hypothesis  to  hold.  This  allows  us  to 
analyze  most  natural  situations,  including  those  listed  above. 

Example  Calculations 

15.  The  previous  description  of  a  variogram  is  admittedly  confusing.  To 
better  explain  the  variogram  function,  a  short  numerical  example  is  forwarded. 
Let  a  linear  arrangement  of  data  be  as  shown  in  Figure  B2 .  At  each  data  loca¬ 
tion,  A  through  G,  a  numerical  value  is  given.  These  locations  and  data 
values  provide  the  information  needed  to  compute  a  variogram. 

16.  This  computation  is  given  in  Table  Bl.  We  begin  by  examining  a 
separation  distance  that  is  close,  but  not  equal,  to  zero.  (For  a  separation 
distance  of  zero,  we  would  be  looking  at  the  difference  between  a  data  value 
with  itself,  which  is  zero.  Hence,  no  useful  information  is  really  obtained 
for  h  equal  to  zero) .  This  example  begins  with  a  separation  distance  of  one 
unit  (h=l)  .  All  data  pairs  separated  by  this  distance  are  obtained  as  shown 
in  Table  Bl;  in  this  case,  six  pairs  are  found.  For  each  pair,  the  difference 
between  values  is  squared.  These  squared  differences  are  then  summed  and 
divided  by  twice  the  number  of  pairs.  This  becomes  the  variogram  value  for  a 
particular  separation  distance. 

17.  Subsequent  to  this,  the  separation  distance  is  incremented  until  a 
distance  is  reached  which  exceeds  the  separation  distance  between  the  farthest 
spaced  pairs.  For  real  data  whose  locations  are  irregularly  spaced,  the  sepa¬ 
ration  distance  h  ,  is  incremented  in  the  same  manner  from  a  small  value  to  a 
larger  value.  An  exception  is  that  h  is  considered  to  be  an  average  value. 
For  example,  for  h  =  2  ,  we  would  find  all  pairs  separated  by  a  distance 
between  1.5  and  2.5  with  an  average  separation  distance  of  2.0. 

18.  The  results  of  a  variogram  calculation  are  displayed  graphically. 

For  the  computation  developed  in  Table  Bl ,  the  graphical  result  is  shown  in 
Figure  B3.  This  variogram  is  seen  to  increase  from  the  origin  as  expected. 
Beyond  a  separation  distance  of  three  units,  however,  the  variogram  begins  to 
decrease.  This  is  occasionally  encountered  with  variogram  computation  and,  in 
this  instance,  shows  that  data  locations  are  similar  at  the  ends  of  this 
linear  array.  It  is  essential  to  remember,  nevertheless,  that  each  set  of 
data  has  a  unique  variogram.  This  is  a  valuable  aspect  of  geostatistics. 
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Table  B1 

Computation  of  a  Variogram  for  the  Data  of  Figure  B1 


Pairs  Difference  (Squared) 


A 

— 

B 

(2.54 

-  2.49)**2 

S 

0.0025 

B 

- 

C 

(2.49 

-  2.43)**2 

= 

0.0036 

C 

- 

D 

(2.43 

-  2. 39)**2 

= 

0.0016 

D 

- 

E 

(2.39 

-  2.41) **2 

= 

0.0004 

E 

- 

F 

(2.41 

-  2.49)**2 

ss 

0.0064 

F_ 

— 

_G 

(2.49 

-  2.56)**2 

* 

0.0049 

N 

ss 

6 

sum 

= 

0.0194 

Y(l)  =  0 

.0194/2(6) 

s 

0.0016 

A 

— 

C 

(2.54 

-  2 . 43) **2 

a 

0.0121 

B 

- 

D 

(2.49 

-  2 . 39)**2 

= 

0.0100 

C 

- 

E 

(2.43 

-  2.41)**2 

0.0004 

D 

- 

F 

(2.39 

-  2 . 49) **2 

= 

0.0100 

E_ 

— 

_G 

(2.41 

-  2.56)**2 

= 

0.0225 

N 

a 

5 

sum 

a 

0.055 

Y(2)  = 

0.055/2(5) 

= 

0.0055 

A 

_ 

D 

(2.54 

-  2 . 39) **2 

0.0225 

B 

- 

E 

(2.49 

-  2.41)**2 

= 

0.0064 

C 

- 

F 

(2.43 

-  2.49)**2 

= 

0.0036 

D_ 

— 

_G 

(2.39 

-  2.56)**2 

* 

0.0289 

N 

= 

4 

sum 

a 

0.0614 

Y (3)  =  0 

.0614/2(4) 

= 

0.0077 

A 

_ 

E 

(2.54 

-  2.41)**2 

= 

0.0169 

B 

- 

F 

(2.49 

-  2.49)**2 

= 

0 

C_ 

— 

_G 

(2.43 

-  2.56) **2 

= 

0.0169 

N 

a 

3 

sum 

a 

0.0338 

Y (4)  =  0 

.0338/2(3) 

= 

0.0056 

A 

aa 

F 

(2.54 

-  2.49)**2 

a 

0.0025 

B_ 

- - 

_G 

(2.49 

-  2 . 56) **2 

= 

0.0049 

N 

r= 

2 

sum 

s 

0.0074 

Y(5)  =  0 

.0074/2(2) 

s 

0.0019 

A 

— 

G 

(2.54 

-  2.56)**2 

a 

0.0004 

N 

= 

1 

Y(6)  =  0 

.0004/2(1) 

s 

0.0002 

No  pairs  found 
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Table  B1  (Concluded) 
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19.  Variogram  estimation  is  fundamental  to  the  understanding  of  spatial 
structure.  Moreover,  because  the  variogram  is  estimated  from  data  represent¬ 
ing  a  particular  spatial  phenomenon,  geostatistics  is  adaptable  to  any  spatial 
situation  without  the  need  for  modification.  The  variogram  is,  as  a  result,  a 
powerful  analytical  procedure. 

20.  The  variogram  is  the  foundation  of  the  regionalized  variables  esti¬ 
mator,  known  as  kriging  (Matheron  1963).  Kriging  is  developed  for  one  primary 
objective:  everywhere  estimates  are  made,  the  variance  of  the  error  of  these 
estimates  is  a  minimum.  Further,  this  noble  objective  must  be  achieved  within 
the  constraint  that  estimation  proceeds  without  changing  the  mean  value  of  the 
spatial  phenomenon.  Hence,  kriging  is  a  constrained,  optimal  estimator. 

21.  A  numerical  example  is  useful  for  showing  how  the  weights,  A^  , 
are  computed  in  the  kriging  system.  Returning  to  Figure  B2,  suppose  an  esti¬ 
mate  is  to  be  made  at  location  H  on  the  basis  of  locations  B  ,  C  ,  and 

D  .  The  geometry  of  this  example  is  simple  and  allows  for  a  clear  description 
of  kriging.  This  does  not  suggest,  however,  that  kriging  always  requires  such 
regular  geometries. 

22.  For  this  example,  a  function  is  required  to  model  the  variogram  of 
Figure  B3.  For  most  applications  of  kriging,  a  spherical  variogram  model  is 
useful.  Such  a  model  has  the  form  of  equation  (B2) . 

23.  In  this  example,  let  equation  (B2)  be  adapted  to  fit  the  variogram 
of  Figure  B3.  For  this  simple  case,  the  following  is  true: 

a.  Nugget  =  0. 

b.  Sill  =  0.0077. 

c.  Range  =  3. 

On  this  basis,  equation  (B2)  becomes 

0.0077 (h/2  -  h3/54)  ,  h  <  3 

y(h)  =  (B8) 

0.0077,  h  £  3 


This  is  the  variogram  model  for  this  example. 

24.  The  procedure  used  to  form  the  kriging  system  for  this  example  is 
shown  in  Table  B2.  Some  explanation  is  required  concerning  the  formation  of 
the  matrix  system.  For  the  matrix,  t°ziZj^  *  last  tow,  (1110)  , 

assures  that  the  sum  of  the  weights  is  equal  to  1.  Further,  with  respect  to 
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Table  B2 

Example  Kriging  Computation 


Step  1.  Form  intersample  covariance  matrix 


Vj 


°7  7  °7  7  °7  7  ^ 

Vl  T2  13 


Z  Z  ®7  Z  ^ 

21  L2LZ  L1L3 


°Z  Z  °7  Z  C7  7  ^ 

j  i  Z3*2  L3  2 


1  0 


where 


^  az  z 

i  i 

sill  -  0.0077 

2)  °7  7  = 

i  j 

a  =  sill  - 

j  i 

T<hu) 

From  this 

i 

j 

hij 

(equation 

1-B 

2-C 

1 

0.0037 

l  =  B 

3-D 

2 

0.0065 

2-C 

3-D 

1 

0.0037 

V3 


0.0040 

0.0012 

0.0040 


Hence 


N 


0.0077 

0.0040 

0.0012 

1 

m 

0.0040 

0.0077 

0.0040 

1 

0.0012 

0.0040 

0.0077 

1 

1 

1 

1 

0 

(Continued) 
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Table  B2  (Continued) 


Step  2.  Form  point  -  sample  covariance  vector 


r 


Z*(H)Z 


A 


Z*(H)Zj 

?Z*(H)Z, 


< 


> 


Z* (H)Z„ 


Z*(H) 

Zi 

Distance,  h 

Y  (h) 

o(h) 

H 

B 

1.414 

0.0050 

0.0027 

H 

C 

1.000 

0.0037 

0.0040 

H 

D 

1.414 

0.0050 

0.0027 

Hence 


°Z*(H)Z 


r  -n 
0.0027 

<  0.0040  ^ 

0.00271 


Step  3.  Solve  the  following  system 


0.0077 

0.0040 

0.0012 

f 

X1 

r  A 

0.0027 

0.0040 

0.0077 

0.0040 

1 

◄ 

»2 

►  < 

0.0040 

0.0012 

0.0040 

0.0077 

1 

»3 

0.0027 

1 

1 

1 

0 

.  1  > 

From  Gauss  Elimination: 


Xj  =  0.289  \2  =  0.422  *3  =  0.289 


V  =  -1.56E-03 


(Continued) 
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Table  B2  (Concluded) 


Step  4 


Step  5 


Compute  the  estimate 

Z*(H)  =  XlZB  +  X2ZC  +  X3Zd 

*  (0.289) (2. 49)  +  (0.422) (2.43)  +  (0.289) (2.39) 

=  2.44 

Compute  the  kriging  variance 

N 

Krig  Var  =  Sill  -  XiaZ*(H)Z  ~  p 

i=l 

=  0.0077  -  [(0.289) (0.0027)  +  (0.422) (0.0040) 
+  (0.289) (0.0027)]  -  (-1.56E-03) 

-  0.0060 

Standard  Deviation  *  0.0060  =  0.078 


(Sheet  3  of  3) 
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this  matrix,  the  unity  values  in  the  fourth  column  are  included  because  the 
Lagrangian  multiplier  helps  relate  the  covariances  of  the  left  hand  side  of 
the  equation  with  those  of  the  right  hand  side.  This  multiplier  assures  the 
minimization  of  the  variance  of  the  error. 

25.  From  this  example,  it  is  observed  that  the  largest  weight  is 
associated  with  location  C  ,  closest  to  the  estimation  location,  H  .  The 
weights,  ,  are  a  function  of  a  separation  distance,  and  closer  points  will 

receive  greater  weight,  in  general,  than  farther  points.  The  estimate 
obtained  in  this  example  is  2.44  with  an  error  having  a  standard  deviation  of 
0.078.  Hence,  kriging  not  only  yields  an  estimate  but  also  provides  an 
estimate  of  its  error. 
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APPENDIX  C:  FRACTAL  DIMENSION  CALCULATIONS 


1.  Table  Cl  presents  the  information  used  to  calculate  the  fractal 
dimensions  of  cross  sections  using  the  method  1  calculation.  Table  C2  pre¬ 
sents  the  information  used  to  calculate  the  fractal  dimensions  of  cross  sec¬ 
tions  using  the  method  2  calculation.  These  tables  are  included  to  show  the 
significant  difference  between  these  techniques  for  the  calculation  of  the 
fractal  dimension  of  rock  surfaces.  The  method  2  calculation  is  the  correct 
procedure . 


Cl 


Table  Cl 

Fractal  Dimension  Summary  (Method  1  Calculation) 


Number 

of  Segment  Counts, 

N 

Fractal 

Cross  Section 

y=305  cm 

y**  100  cm 

y=20  cm 

Dimension 

AB 

4 

13 

65 

1.023 

AF 

4 

12 

62 

1.006 

AG 

5 

15 

78 

1.008 

AH 

4 

12 

64 

1.018 

AJ 

5 

15 

79 

1.013 

AN 

3 

9 

46 

1.002 

AP 

3 

10 

50 

1.033 

AQ 

3 

10 

51 

1.040 

AR 

4 

12 

63 

1.012 

AS 

3 

10 

51 

1.040 

AT* 

2 

7 

36 

1.061 

AU* 

1 

5 

25 

1.181 

AW* 

2 

8 

41 

1.109 

AX* 

2 

6 

33 

1.029 

AY* 

2 

7 

39 

1.090 

AZ* 

1 

4 

24 

1.166 

BB 

8 

24 

123 

1.003 

BC 

8 

24 

123 

1.003 

BD 

8 

24 

125 

1.009 

BE 

7 

23 

119 

1.040 

BG 

6 

19 

100 

1.033 

BH 

7 

22 

113 

1.021 

BI 

6 

19 

96 

1.018 

BJ 

6 

20 

101 

1.036 

*  The  cross-sectional  length  of  sections  AT  through  AZ  is  less  than  that  of 
sections  AB  through  AS.  Hence,  the  segment  length,  y  *  305  cm  ,  was  too 
large.  As  a  result,  the  use  of  whole  numbers  to  represent  the  number  of 
segment  counts,  N  ,  was  incorrect. 


Table  C2 

Fractal  Dimension  Summary  (Method  2  Calculation) 


Cross  Section 

Number 
y=305  cm 

of  Segment  Counts, 
y=100  cm 

N 

y*20  cm 

Fractal 

Dimension 

AG 

4.951 

15.200 

77.000 

1.0072 

AH 

4.164 

12.750 

65.600 

1.0119 

AJ 

4.974 

15.480 

77.100 

1.0060 

AN 

3.013 

9.200 

46.100 

1.0012 

AP 

3.328 

10.170 

50.200 

1.0009 

AQ 

3.302 

10.160 

51.200 

1.0061 

AR 

4.000 

12.490 

62.800 

1.0107 

AT 

2.285 

7.140 

36.800 

1.0200 

AW 

2.551 

8.140 

41.850 

1.0268 

AX 

2.049 

6.290 

32.900 

1.0189 

C3 


APPENDIX  D:  NUMERICAL  VARIOGRAM  RESULTS 


1.  This  appendix  presents  the  numerical  results  used  to  plot  the 
variograms  which  were  presented  in  Part  III.  These  results  are  presented  to 
verify  each  variogram  plot  and  to  provide  numerical  precision  beyond  what  can 
be  interpreted  from  each  plot.  Only  results  for  the  variograms  of  Part  III 
are  presented. 
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Table  D1 

Variogram  Results  for  Cross  Section  BG 


h(cm) 

Y(h) 

No.  of 
Pairs  (1 

0  - 

50 

6.00 

48 

51  - 

100 

16.96 

49 

101  - 

150 

30.32 

53 

151  - 

200 

49.38 

42 

201  - 

250 

73.94 

49 

251  - 

300 

117.47 

44 

301  - 

350 

164.15 

40 

351  - 

400 

237.89 

41 

401  - 

450 

311.95 

39 

451  - 

500 

377.23 

28 

501  - 

550 

416.96 

38 

551  - 

600 

517.00 

41 

601  - 

650 

671.00 

26 

651  - 

700 

722.03 

19 

701  - 

750 

821.73 

22 

751  - 

800 

915.75 

20 

801  - 

850 

1057.55 

22 

851  - 

900 

1161.12 

17 

901  - 

950 

1321.39 

18 

951  - 

1000 

1474.82 

17 

Sample  Variance  *  427.81 


Table  D2 

Variogram  Results  for  Cross  Section  AF 


No.  of 


h(cm) 

Y(h) 

Pairs  (N) 

0 

- 

80 

85.71 

28 

81 

- 

160 

232.26 

31 

161 

- 

240 

242.19 

32 

241 

- 

320 

258.33 

24 

321 

- 

400 

364.58 

24 

401 

- 

480 

707.14 

21 

481 

- 

560 

816.67 

21 

561 

- 

640 

871.88 

16 

641 

- 

720 

1125.00 

18 

721 

- 

800 

1635.00 

20 

801 

- 

880 

2261.77 

17 

881 

- 

960 

2335.71 

14 

961 

- 

1040 

2512.50 

16 

1041 

- 

1120 

2875.00 

16 

1121 

- 

1200 

3290.00 

15 

1201 

- 

1280 

4783.33 

18 

1281 

- 

1360 

5525.00 

10 

1361 

- 

1440 

4075.00 

4 

1441 

— 

1520 

6658.33 

6 

Sample  Variance  =  1525.64 
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Table  D3 

Variogram  Results  for  Cross  Section  AL 


h(cm) 

Y(h) 

No.  of 
Pairs  (N) 

0-40 

67.38 

71 

41  -  80 

207.85 

80 

81  -  120 

418.24 

76 

121  -  160 

600.19 

75 

161  -  200 

848.87 

76 

201  -  240 

912.32 

74 

241  -  280 

1165.48 

65 

281  -  320 

1217.34 

64 

321  -  360 

1663.28 

58 

361  -  400 

1679.67 

61 

401  -  440 

1581.22 

49 

441  -  480 

1279.12 

57 

481  -  520 

1087.29 

42 

521  -  560 

815.43 

42 

561  -  600 

684.22 

37 

601  -  640 

1034.19 

31 

641  -  680 

1587.79 

28 

681  -  720 

2254.17 

24 

721  -  760 

2758.40 

30 

761  -  800 

4054.22 

18 

Sample  Variance  =  1392.00 
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Table  D4 

Variogram  Results  for  Cross  Section  AS' 


No.  of 

h(cm)  y(h)  Pairs  (N) 


0 

- 

40 

15.57 

54 

41 

- 

80 

72.24 

59 

81 

- 

120 

137.14 

61 

121 

- 

160 

204.16 

55 

161 

- 

200 

279.22 

48 

201 

- 

240 

401.56 

47 

241 

- 

280 

411.65 

44 

281 

- 

320 

548.44 

48 

321 

- 

360 

775.49 

41 

361 

- 

400 

944.30 

33 

401 

- 

440 

870.51 

45 

441 

- 

480 

1029.90 

34 

481 

- 

520 

1030.64 

36 

521 

- 

560 

978.05 

28 

561 

- 

600 

885.94 

26 

601 

- 

640 

912.30 

23 

641 

- 

680 

875.28 

16 

681 

- 

720 

1083.47 

17 

721 

- 

760 

1511.42 

13 

761 

— 

800 

2043.61 

9 

Sample  Variance  =  649.39 
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Table  D5 

Variogram  Results  for  Cross  Section  BA 


h(cm) 

Y(h) 

No.  of 
Pairs  (N) 

0  - 

50 

26.60 

37 

51  - 

100 

93.65 

40 

101  - 

150 

129.00 

34 

151  - 

200 

127.79 

29 

201  - 

250 

108.22 

36 

251  - 

300 

75.55 

31 

301  - 

350 

41.86 

29 

351  - 

400 

49.67 

18 

401  - 

450 

78.40 

20 

451  - 

500 

45.89 

19 

501  - 

550 

39.69 

13 

551  - 

600 

44.00 

13 

601  - 

650 

55.71 

14 

651  - 

700 

128.43 

14 

701  - 

750 

86.00 

16 

751  - 

800 

32.77 

13 

801  - 

850 

60.55 

11 

851  - 

900 

28.00 

11 

901  - 

950 

20.00 

9 

951  - 

1000 

103.14 

7 

Sample  Variance  *  105.47 
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Table  D6 

Variogram  Results  for  Cross  Section  BB 


No.  of 


h(cm) 

Y(h) 

Pairs 

0  - 

50 

10.52 

84 

51  - 

100 

28.45 

71 

101  - 

150 

38.96 

71 

151  - 

200 

42.11 

68 

201  - 

250 

49.86 

69 

251  - 

300 

66.79 

58 

301  - 

350 

80.38 

63 

351  - 

400 

85.70 

54 

401  - 

450 

91.88 

69 

451  - 

500 

110.67 

63 

501  - 

550 

175.91 

46 

551  - 

600 

196.44 

45 

601  - 

650 

256.33 

55 

651  - 

700 

292.09 

47 

701  - 

750 

323.66 

41 

751  - 

800 

394.04 

56 

801  - 

850 

436.00 

50 

851  - 

900 

579.28 

39 

901  - 

950 

509.30 

37 

951  - 

1000 

512.12 

34 

Sample  Variance  =  263.77 
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Table  D7 

Variogram  Results,  Cross  Section  From  Patton  and  Deere  (1970) 


No.  of 

h(cm)  y (h)  Pairs  (N) 


0 

- 

8 

0.25 

269 

9 

- 

16 

0.75 

342 

17 

- 

24 

1.68 

326 

25 

- 

32 

2.53 

310 

33 

- 

40 

3.35 

294 

41 

- 

48 

3.76 

274 

49 

- 

56 

3.84 

258 

57 

- 

64 

3.71 

242 

65 

- 

72 

3.28 

226 

73 

- 

80 

2.63 

210 

81 

- 

88 

1.90 

194 

89 

- 

96 

1.14 

178 

97 

- 

104 

0.78 

162 

105 

- 

112 

0.86 

146 

113 

- 

120 

1.70 

130 

121 

- 

128 

3.04 

114 

129 

- 

136 

4.32 

98 

137 

- 

144 

5.11 

79 

145 

- 

152 

5.17 

62 

153 

- 

160 

4.13 

55 

Sample  Variance  =  2.42 
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Table  D8 

Variogram  Results  for  Cross  Section  BG  (Filtered) 


No.  of 


_h_ 

y(h) 

Pairs  (N) 

0 

- 

20 

1.15 

15 

21 

- 

40 

2.65 

23 

41 

- 

60 

2.94 

20 

61 

- 

80 

6.29 

21 

81 

- 

100 

6.62 

18 

101 

- 

120 

5.45 

20 

121 

- 

140 

7.47 

20 

141 

- 

160 

8.87 

19 

161 

- 

180 

9.11 

19 

181 

- 

200 

6.32 

17 

201 

- 

220 

4.92 

21 

221 

- 

240 

4.54 

19 

241 

- 

260 

4.17 

19 

261 

- 

280 

4.96 

18 

281 

- 

300 

3.92 

16 

301 

- 

320 

2.66 

14 

321 

- 

340 

2.16 

18 

341 

- 

360 

4.16 

16 

361 

- 

380 

7.55 

15 

381 

— 

400 

8.16 

22 

Sample  Variance  =  5.92 
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APPENDIX  E:  FORTRAN  PROGRAM  FOR  CALCULATION 
OF  THE  FRACTAL  DIMENSION 


1.  This  appendix  presents  a  FORTRAN  program  for  the  calculation  of  the 
fractal  dimension  of  a  one-dimensional  irregular  line.  This  program  was 
developed  to  analyze  the  string  line  cross  sections  presented  in  Part  V  of 
this  report.  A  short  user's  guide  follows. 


User's  Guide 


2.  Input  to  this  program  is  afforded  through  a  free  format  structure; 
this  facilitates  data  entry.  Input  consists  of  three  sections: 

a.  Record  1:  READ (5,*)  K 

where  K  is  the  number  of  step  sizes,  y  ,  to  be 
used  to  calculate  the  fractal  dimension.  For 
Table  6,  K  is  4. 

b.  Record  2:  READ(5 , *)  (SL(I),  I  -  1,K) 

where  SL  is  a  vector  which  contains  the  step  sizes, 
y  .  For  Table  6,  these  values  are  0.5,  1.0,  1.5, 
and  2.0. 

c.  Record (s)  3:  Data  entry;  one  record  is  required  for  each 

coordinate  location.  For  a  50-ft  cross  section, 
for  example,  sampled  every  0.5  ft,  100  records 
are  required.  Each  record  is  of  the  form: 

READ (5,*)  ELEVATION,  X-COORDINATE 
The  last  record  is 

ELEVATION  =  -999.0 
X-COORDINATE  =  -999.0 


El 


'O  00  -4  ON  Xr  Ui  ro 


PROGRAM  FRACTAL  FORTRAN 


C 

C 

C 

C 


c 

«* 

»* 

c 

*« 

»« 

c 

•« 

«# 

c 

«« 

«« 

c 

»* 

## 

c 

»* 

«# 

c 

«* 

«« 

c 

«« 

THIS 

PROGRAM 

PERFORMS  TO  DETERMINE  THE  PARTICULAR 

«* 

c 

«« 

NUMBER 

OF  STEPS,  EN,  USED  TO  APPROXIMATE  THE 

c 

#* 

GLOBAL 

LENGTH 

FOR  FRACTAL  DIMENSION  CALCULATION. 

«# 

c 

»» 

AND 

TO 

CALCULATE  THE  FRACTAL  DIMENSION. 

#* 

c 

## 

«# 

c 

#« 

GUIDE  TO  DATA  INPUT 

## 

c 

## 

## 

c 

«* 

CARD 

1 

K, 

NUMBER  OF  PARTICULAR  SEGMENTS 

#* 

c 

«« 

c 

## 

CARD 

2 

SL, 

LENGTH  VALUES  OF  SEGMENTS 

c 

#« 

c 

«# 

CARD 

3 

DATA  ENTRY 

#4f 

c 

#* 

Y 

1ST  PART,  ELEVATION 

c 

## 

X 

2ND  PART,  X-COORDINATE 

## 

c 

## 

c 

LAST 

CARD 

2ND  PART  (X)  =  NEGATIVE  VALUE 

c 

## 

C  ###*#####***####*###*######*#**#**#* «**«««**«****#«»«**««««* 

c 

1  IMPLICIT  REAL*8  (A-H.O-Z) 

COMMON/COORO/  X(1000),  Y(TOOO) 

COMMON/O I  ST/  DT(IOOO),  SL( 10),  EN(10) 

COMMON/CONST/  D(10,10),  XX(10) 

COMMON/CENTER/  XC.  YC,  XI,  Y1 

DIMENSION  TEMP( 10,11),  F { 10),  CK(10),  LOC(IO),  B(10) 

C 

DO  5  J  =  1,  10 
DO  5  1=1,  10 
D(  l,J)  =  0.0 

10  XX ( I )  =  0.0 

11  5  CONTINUE 

12  DO  10  I  =  1,  1000 

13  X( I )  =  0.0 

14  Y(  I  )  =  0.0 

15  DT ( I )  =  0.0 

16  10  CONTINUE 

C 

C. . . .  READ  THE  OATA 
C 


E2 


17  READ(5,*)  K 

18  READ( 5, # )  ( SL( I ) ,  I  =  1,  K) 

19  KOUNT  =  0 

20  20  READ( 5,  * )  DUM1 ,  0UM2 

21  I F( 0UM2 . LT . 0 . 0 )  GO  TO  30 

22  KOUNT  =  KOUNT  +  1 

23  Y( KOUNT)  =  DUM1 

24  X( KOUNT)  =  DUM2 

25  GO  TO  20 

C 

C _  CALCULATE  ALL  DISTANCES  BETWEEN  TWO  POINTS 

C  AND  FIRST  CENTER  OF  CIRCLE  (XC,  YC) 

C 

26  30  ND  =  KOUNT 

27  N  =  ND  -  1 

28  DO  40  I  =  1,  N 

29  XD  =  X(  I )  -  X( 1+1  ) 

30  YD  =  Y( I )  -  Y( 1+1 ) 

31  DT ( I )  =  DSQRT ( XO  #  XD  +  YD  *  YD) 

32  40  CONTINUE 

C 

C -  LOOP  OF  NUMBERS  OF  DIFFERENT  SEGMENTS  (R) 

C 

33  DO  10000  I  I  =  1,  K 

34  SUM  =0.0 

35  NKOUNT  =  0 

36  R  =  SL(  I  I  ) 

37  WR I TE( 6, 50  )  R 

38  50  FORMAT( 1H1 ,  5(/),  10X,  '*»**###*  SEGMENT  LENGTH  =  F5.2, 

&  '  . **'.  3(  / )  ) 

39  WR I TE( 6, 60 )  ND 

40  60  FORMAT ( 6X,  '  THE  NUMBER  OF  DATA  =  ’ ,  14,  4(/), 

&  1 9X,  'INTERCEPT  POINTS',  /,  19X,  '--- . . ' , 

8c  //,  13X,  'NUMBER',  4X,  'X-COORD',  4X,  'Y-COORD',  // ) 

41  XC  =  X( 1 ) 

42  YC  =  Y( 1 ) 

C 

C....  FIND  DT( M )  USING  COMPARISION  OF  THE  SUM  OF  EACH  DISTANCE 
C 

43  J  =  1 

44  1000  DO  70  I  =  J,  N 

45  SUM  =  SUM  +  DT( I ) 

46  I F( SUM.GE.R)  GO  TO  80 

47  I F ( I . GE . N )  GO  TO  3000 

48  70  CONTINUE 

49  80  M  =  I 

C 

C _  CALL  SUBROUTINE  TO  CALCULATE  THE  INTERCEPT  POINT 

C 

50  2000  CALL  CROSS  (M,  R) 
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51 

52 

53 

54 

55 

56 

57 

58 

59 


60 

61 

62 

63 

64 

65 

66 


67 

68 
69 


70 

71 

72 

73 

74 

75 

76 

77 


78 

79 


C 

C. 

C 

C 


90 


100 

C 

C.  . 

C 

C 


C 

C.  .  . 

C 

C 

3000 


CHECK  WHETHER  THE  INTERCEPT  POINTS  ARE  INTERVAL  BETWEEN  TWO 
POINTS,  (X(M),  Y( M) )  &  (X(M+1),  Y(M+1)) 

I F(X1 .GE.X(M)  .AND.  XI . LE . X( M+1 ) )  XC  =  XI 
YY  =  (Y(M)  -  Y1)  *  ( Y( M+1 )  -  Y1 ) 

IF(YY.LE.O.O)  GO  TO  90 
M  =  M  +  1 
GO  TO  2000 
YC  =  Y1 

NKOUNT  =  NKOUNT  +  1 

WR I TE( 6, 100 )  NKOUNT,  XC,  YC 

FORMAT ( 1 3X,  14,  4X,  F8.3,  3X,  F8.3) 

CALCULATE  THE  DISTANCE  (ON)  BETWEEN  THE  INTERCEPT  POINT  (NEW 
CENTER  OF  CIRCLE)  AND  NEXT  POINT  (X(M+1),  Y(M+1)) 

XN  =  XC  -  X( M+1 ) 

YN  =  YC  -  Y( M+1 ) 

DN  =  DSQRT ( XN  *  XN  +  YN  *  YN ) 

SUM  =  DN 

IF(SUM.GE.R)  GO  TO  2000 
J  =  M  +  1 
GO  TO  1000 


CALCULATE  THE  REMAINING  DISTANCE  (DR) 

FROM  CENTER  OF  CIRCLE  AT  LAST  PORTION 

XR  =  XC  -  X(ND) 

YR  =  YC  -  Y(ND) 

DR  =  OSQRT(XR  *  XR  +  YR  *  YR ) 


CALCULATE  THE  EXACT  NUMBER  OF  STEPS  (EN) 


EN(  I  I )  =  NKOUNT  +  DR  /  R 
WRITE(6, 110)  EN(  I  I  ) 

110  FORMAT ( 5( / ) ,  5X,  'NUMBER  OF  STEPS  (N)  =  ',  F10.4) 

10000  CONTINUE 

WRITE(6, 120) 

120  FORMAT ( 1  HI ,  5(/),  10X,  'SEGMENT  LENGTH',  4X,  ' NUMBER  OF  STEPS’, 

&  /,  10X,  ( . ',  4X,  ' . . — ',  //) 

WR I TE( 6, 130)  ( SL( I ) ,  EN(I),  I  =  1,  K) 

130  FORMAT ( 1 5X,  F5.L,  11X,  F8.4,  /) 

C 

C _  CALCULATION  OF  FRACTAL  DIMENSION 

C  USING  SEGMENT  LENGTH  AND  NUMBER  OF  STEPS 

C -  CALCULATION  OF  LOG  10  SCALE  FOR  DATA 

C 


DO  5100  I  =  1,  K 

SL( I )  =  DL0G10  (SL( I )) 
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80 

EN( 1 )  =  DL0G10  ( EN(  1  )  ) 

81 

5100 

CONTINUE 

C 

C.  .  .  . 

FORM  THE  MATRIX  jDj 

C 

82 

D(1,1)  =  K 

83 

DO  5200  1  =  1,  K 

84 

D( 1,2)  =  D( 1 , 2  )  +  SL( 1 ) 

85 

0(2,1)  =  0(1,2) 

86 

0(2,2)  =  D( 2, 2 )  ♦  SL(  1 ) 

*  SL( 1 ) 

87 

5200 

CONTINUE 

C 

C.  .  .  . 

FORM  THE  VECTOR  )XXj 

C 

88 

DO  5300  1  =  1  ,  K 

89 

XX(  1 )  =  XX( 1  )  +  EN(  1  ) 

90 

XX(  2  )  =  XX(  2  )  +  SL(  1  )  * 

EN(  1  ) 

91 

5300 

CONTINUE 

C 

C.  .  .  . 

SOLVE  THE  LINER  COEFFICIENTS,  VECTOR 

C 

92 

CALL  EQSOLD  (2,  F) 

C 

93 

FD  =  -  F (  2 ) 

94 

WR 1 TE( 6, 5400  )  FD 

95 

5400 

FORMAT ( 5 ( / ) ,  5X,  ‘FRACTAL 

DIMENSION  = 

96 

STOP 

97 

END 

c 

C _  SUBROUTINE  TO  CALCULATE  THE  INTERCEPT  POINTS 

C  BETWEEN  STRAICHT  LINE  AND  CIRCLE. 

C 


5(/)) 


98 

99 
100 
101 
102 


103 

104 

105 

106 


SUBROUTINE  CROSS  (M,  R) 

IMPLICIT  REAL*8  (A-H.O-Z) 

COMMON/COORD/  X(1000),  Y(1000) 

COMMON/D  I  ST /  DT(IOOO),  SL(10),  EN(10) 
COMMON/CENTER/  XC,  YC,  XI,  Y1 
C 

C _  CALCULATE  THE  COEFFICIENTS  OF  STRAIGHT  LINE 

C  USING  (  X(M),  Y(M)  )  &  (  X(M+1),  Y(M-H)  ) 

C  Y  =  YM  /  XM  *  (X  -  X( M ) )  ♦  Y(M) 

C 

XM  =  X( M+1 )  -  X(M) 

YM  =  Y(M+1 )  -  Y(M) 

A1  =  YM  /  XM 

A2  =  ( — A 1 )  *  X(M)  +  Y(M) 

C 

C _  CALCULATE  THE  COEFFICIENTS  OF  2ND-0RDER  EQUATION 

C  USING  STRAIGHT  LINE  AND  CIRCLE 

C 
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C  CIRCLE  :  (X  -  XC)**2  +  (Y  -  YC)**2  =  R**2 

C 

C  2N0-0R0ER  EQ  :  A  *  X**2  +  B*X  +  C  =  0 

C 


107 

A  =  A1  *  A1  +  1 .0 

108 

B  =  2.0  *  A1  *A2  -  2.0  *  XC  -  2.0  *A1  *  YC 

109 

CC  =  XC  *  XC  +  YC  *  YC  -  R  *  R 

110 

Z' 

C  =  A2  *  A2  -  2.0  *  A2  *  YC  +  CC 

Vs 

c. . . 

.  CALCULATE  THE  SOLUTION  OF  2ND-ORDER  EQUATION 

c 

n 

USING  X  =  ( -B  +0R-  SQRT ( B*B-4#A*C ) )  /  (2  *  A) 

1 1 1 

Vs 

OR  =  B  *  B  -  4.0  *  A  *  C 

112 

DD  =  DSQRT(DR) 

113 

XI  =  ( -B  +  DD)  /  (2.0  *  A) 

114 

VI  =  A1  *  XI  +  A2 

115 

RETURN 

116 

r 

END 

117 

Vs 

p 

SUBROUTINE  EQSOLD  <N,  F) 

Vs 

c. .  . 

.  SUBROUTINE  TO  PERFORM  EQUATION  SOLUTION 

c 

r> 

USING  GAUSS  EL  1 M 1 NA 1  ON 

118 

Vs 

IMPLICIT  REAL*8  (A-H.O-Z) 

119 

COMMON/CONST/  0(10,10),  XX(10) 

120 

DIMENSION  TEMP( 10,11),  F ( 10),  CK{10),  LOC(IO),  B(10) 

121 

MP  =  N  +  1 

122 

DO  20  1  =  1,  N 

123 

DO  10  J  =  1,  N 

124 

TEMP( 1 , J )  =  D{ 1 , J ) 

125 

10 

CONTINUE 

126 

TEMP(  1 ,MP)  =  XX( 1  ) 

127 

20 

CONTINUE 

128 

DO  30  1  =  1,  N 

129 

CK( 1 )  =  0.0 

130 

30 

CONTINUE 

131 

DO  100  1  =  1,  N 

132 

p 

1  P  =  1  +  1 

c. . . 

p 

FIND  MAX.  NUMBER  IN  THE  ITH  COLUMN 

133 

Vs 

AMAX  =  0.0 

134 

DO  40  K  =  1,  N 

135 

1 F ( AMAX  -  DABS( TEMP( K, 1 ) ) )  35,  40,  40 

136 

35 

1 F ( CK( K )  )  36,  36,  40 

137 

36 

L0C( 1 )  =  K 

138 

AMAX  =  DABS( TEMP( K,  1  )  ) 

139 

40 

CONTINUE 

140 

1 F ( AMAX  .EQ.  0.0)  GO  TO  100 
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c 

C _  MAX.  ELEMENT  IN  THE  ITH  COLUMN  IS  A( L, I ) 

C 


141 

L  =  LOC( 1  ) 

142 

CK( L)  =  1.0 

143 

DO  50  J  =  1,  N 

144 

IF(L-J)  41,  50,  41 

145 

41 

G  =  -TEMP(J, 1 )  /  TEMP( L, 1 ) 

146 

DO  45  K  =  IP,  MP 

147 

45 

TEMP( J, K)  =  TEMP( J , K )  +  G  *  TEMP(L,K) 

148 

50 

CONTINUE 

149 

100 

CONTINUE 

150 

DO  200  1  =  1,  N 

151 

L  =  LOC(  1  ) 

152 

F(  1  )  =  TEMP(L.MP)  /  TEMP(  L,  1  ) 

153 

200 

CONTINUE 

154 

999 

CONTINUE 

155 

RETURN 

156 

END 
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APPENDIX  F:  VARIOGRAMS  FOR  STRING  LINE  CROSS  SECTIONS 


1.  This  appendix  presents  the  variograms  for  the  string  line  cross 
sections  introduced  in  Part  V  of  this  report.  Because  these  variograms  are  of 
minor  significance  to  the  conclusions  of  this  report,  no  further  discussion  of 
them  is  made  other  than  their  presentation  in  this  appendix. 
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Figure  F 4.  Variogram  for  the  Island  Notch  joint  surface 
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Figure  F5.  Variogram  for  the  DS  + 
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Figure  F8.  Variogram  for  the  Dunn 
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Figure  F9.  Variogram  for 
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Figure  F10.  Variogram  for  the  Island  Notch  joint/bedding  plane 
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Figure  Fll.  Variogram  for  the  Island  Notch  bedding  surface 
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Figure  FI  3.  Variogram  for  the  Dunn  Cree 
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Figure  F14.  Variogram  for  the  Old 
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Figure  F15.  Variogram  for  the  Backus 


